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SECTION  I 


INTRODUCTION 

The  VIBRA-4  digital  computer  program  (Ref.  1)  analyzes  the  gross 
structural  response  (bending  moments  and  accelerations)  of  an  aircraft 
exposed  to  a  nuclear  blast.  For  aircraft  flying  at  supersonic  speeds, 
the  aerodynamic  loads  are  determined  in  VIBRA-4  by  a  lifting  surface 
method  that  utilizes  the  indicial  aerodynamic  influence  coefficient  (AIC) 
approach  which  is  an  extension  of  the  well-known  Mach  box  method  that 
has  been  used  in  supersonic  flutter  analyses  for  many  years.  For  sub¬ 
sonic  speeds,  however,  VIBRA-4  has  used  an  indicial  strip  theory  based 
on  two-dimensional  theory  or  empirical  results  derived  from  shock-tube 
measurements.  The  subsonic  Doublet-Lattice  Method  (DLM,  Refs.  2  and  3) 
has  received  wide  acceptance  since  its  appearance  in  1968  for  subsonic 
flutter  analysis  because  of  its  generality  regarding  configurations. 
Current  extensions  (Refs.  3,  4  and  5)  account  for  wing-body  interference 


1.  Hobbs,  N.P.,  Zartarian,  G.,  and  Walsh,  J.P.,  A  digital  Computer  Pro¬ 
gram  for  Calculating  the  Bla*t  Pe*pon*e  of  A ircraft  to  Nuclear  Exploiioni , 

Air  Force  Weapons  Laboratory,  Report  No.  AFWL-TR-70-140,  Vol .  I,  April  1971. 

2.  Albano,  E.,  and  Rodden,  W.P.,  "A  doublet- Lattice  Method  for  Calculating 
Li.ft  diAtributicn*  on  Oiciilatinq  Surface*  in  Subsonic  Flow*,"  AIAA  J.,  Vol. 
7,  No.  2,  pp.  279-285,  and  No.  11,  p.  2192,  1969. 

3.  Giesing,  J.P.,  Kalman,  T.P.,  and  Rodden,  W.P.,  Subionic  Unsteady  Amp- 
dynamic*  for  General  Configuration*;  Part  L  -  direct  Application  of  the 
Nonplanar  doublet- Lattice  Method,  Air  Force  Flight  Dynamics  Labortory, 

Report  No.  AFFDL-TR-71 -5,  Part  I,  November  1971. 

4.  Giesing,  J.P.,  Kalman,  T.P.,  and  Rodden,  W.P.,  Sub-ionic  Unsteady  Aero¬ 
dynamic*  for  General  Con  figuration*;  Part  II  -  Application  of  the  doublet- 
Lattice  Method  and  the  Method  of  Image*  to  Lifting -Surf ace/Body  Interference , 

Air  Force  Flight  Dynamics  Laboratory,  Report  No.  AFFDL-TR-71 -5,  Part  II, 

April  1972. 

5.  Giesing,  J.P.,  Kalman,  T.P.,  and  Rodden,  W.P.  ,'‘Sub*onic  Steady  and 
Oicillat.ory  Aerodynamic*  for  Multiple  Interfering  Wing*  and  Bo  die* , "  J . 
Aircraft,  Vol.  9,  No.  10,  pp.  693-702,  1972 
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and  permit  the  determination  of  the  aerodynamic  loading  on  an  entire 
aircraft  in  oscillatory  motion  or  in  a  harmonic  gust  field  (Ref.  6). 

The  generality  of  the  DLM  makes  it  a  desirable  replacement  for  the  present 
subsonic  capability  of  VIBRA-4.  The  calculation  of  transient  airloads 
from  an  oscillatory  theory  can  only  be  made  by  the  use  of  Fourier  trans¬ 
form  methods;  however,  the  feasibility  of  calculating  the  transient  blast 
loads  by  the  DLM  and  Fourier  techniques  has  already  been  demonstrated  by 
Zartarian  (Ref.  7),  although  it  is  necessary  to  add  Piston  Theory  (Ref.  8) 
to  supplement  the  DLM  calculations  at  high  reduced  frequency  as  was  pointed 
out  by  Zartarian. 

The  VIBRA-6  computer  program,  documented  here,  is  an  extension  of 
earlier  versions  of  VIBRA.  The  purpose  of  this  report  is  to  document  the 
modifications  to  adapt  the  Doublet-Lattice  Method  and  Piston  Theory  to  the 
requirements  of  the  new  VI3P.A-5  computer  program.  Only  the  modifications 
are  described  here  since  the  documentation  for  the  DLM  and  Method  of  Images 
and  for  Piston  Theory  is  quite  extensive.  A  brief  summary,  however,  is 
presented  in  Volume  I  of  this  report.  The  modifications  include  addition 
of  a  travelling  gust  field,  improvements  in  the  aerodynamic  influence  coef¬ 
ficients,  a  compressible  Slender  Body  Theory,  and  changes  in  the  aerodynamic 
load  output.  This  Volume  also  serves  as  a  User's  Manual  since  sections  are 
included  that  1)  discuss  considerations  for  modelling  the  entire  aircraft, 

2)  provide  the  program  input  and  3)  present  an  example  problem.  This  Volume 
also  presents  a  description  of  subroutines  used  in  the  aerodynamic  module. 


6.  Giesing,  J.P.,  Rodden.W.P.,  and  Stahl,  B.,  "SeanA  F unction  and  Lifting 
Surface  Theory  fan  Hanmoni c  GuAt  PieXdA"  ,J.  Aircraft,  Vol  No.  3,  pp.  252- 
255,  1970. 

7.  Zartarian,  G. ,  Application  of  the  doublet- Lattice.  Method  fan  VetenmXna- 
tion  of  ElaAt  Loads  on  Lifting  SuA faces  aX  Subsonic  Speeds,  Air  Force 
Weapons  Laboratory,  Report  No.  AFWL-TR-72-207,  January  1973. 

8.  Ashley,  H.,  and  Zartarian,  G.,  "Pi  A  ton  The.oA.ij  -  A  new  A  eAodijnom-ic  Toot 
fen  the  Acau v las ti clan" , J.  Aero.Sci .  ,Vol .  23, No.  12,  pp.  1109-1118,  1956 
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SECTION  II 


THEORETICAL  DEVELOPMENT 

1.  DOWNWASH  BOUNDARY  CONDITIONS 

The  loads  imposed  on  a  vehicle  flying  through  an  atmospheric  gust  or 
a  blast  wave  arise  from  three  sources.  The  first  source  is  the  aerodynamic 
angle  of  attack  induced  by  the  gust  or  blast  wave,  the  second  source  is 
the  angle  of  attack  induced  by  the  response  motion,  and  the  third  is  the 
impulsive  diffraction  loads  induced  as  the  pressure  wave  sweeps  across 
the  structure.  It  is  generally  assumed  that  the  gust  field  or  blast  wave 
is  known  and  is  independent  of  the  response,  and  that  the  additional 
aerodynamic  loads  induced  by  the  motion  must  be  determined  by  solution 
of  the  aeroelastic  equations  of  motion.  The  aerodynamic  loads  are  not 
known  in  general  for  arbitrary  time-dependent  downwash  distributions 
(from  motion  or  the  gust)  but  only  for  harmonic  time  dependence  of  the 
downwash.  The  response  problem  is  therefore  easily  formulated  in  the 
frequency  domain  and  the  time  history  of  response  follows  from  an  inverse 
Fourier  transformation.  The  specific  intermediate  problems  are:  1)  the 
determination  of  the  harmonic  downwash  caused  by  a  blast  wave  from  an 
arbitrary  burst  orientation,  and  2)  the  harmonic  downwash  caused  by  motion 
of  the  rigid  body  or  vibration  modes  of  the  vehicle. 

2.  DOWNWASH  INDUCED  BY  A  TRAVELLING  HARMONIC  GUST 

An  important  feature  of  the  travelling  gust  is  that  the  aerodynamic 
loads  do  not  depend  simply  on  the  relative  velocity  between  the  gust  and 
the  vehicle.  The  downwash  is,  however,  determined  by  the  relative  velocity 
but  the  basic  analysis  (Doublet-Lattice  formulation)  depends  only  on  the 
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velocity  U  of  the  vehicle  through  the  atmosphere  (assumed  to  be  fixed). 

A  blast  wave  is  a  special  case  of  a  travelling  gust.  A  travelling  gust 
is  usually  considered  to  have  a  transverse  velocity  profile  that  is 
perpendicular  to  the  direction  of  travel.  In  the  analysis  to  follow,  the 
gust  velocity  G  is  the  material  velocity  of  the  fluid  and  is  in  the 
direction  of  travel  of  the  blast  wave,  and  the  velocity  Vg  is  the  velocity 
of  propagation  of  the  pressure  wave. 


For  introductory  purposes  we  consider  the  simplest  case  of  a  vertical 
blast  wave  encountering  a  horizontally  moving  aircraft.  We  denote  the 
gust  velocity  by  G(0  where  t  denotes  the  penetration  distance  into  the 
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gust;  i.e.,  the  distance  normal  to  the  gust  wave  front.  Figure  1  shows 
the  initial  positions  of  the  gust  field  and  the  aircraft  and  the  positions 
at  a  time  t  after  penetration.  From  Figure  1  we  find  the  penetration 
distance  is 

C  -  Vgt  -  (z  +  za)  (1) 

where  za  is  the  coordinate  of  the  aircraft  reference  point  in  the  global 
xa  -  za  coordinate  system,  and  z  is  the  coordinate  of  the  point  of  interest 
relative  to  the  reference  point  in  the  local  x-z  coordinate  system. 

In  this  elementary  case,  the  penetration  distance  is  independent  of 
the  aircraft  velocity.  As  discussed  in  the  Introduction,  the  oscillatory 
aerodynamic  lifting  surface  analysis  requires  a  harmonic  downwash  ampli¬ 
tude  which  may  be  determined  from  a  harmonic  gust  profile  with  wave  length 


x.  Letting 

G(0  =  G0  cos  2n/x  (2a) 

=  Re  [w(x,y ,z)  exp(iwt)]  (2b) 

we  find  the  frequency  is 

oj  =  2tt  V  q/ X  (3) 

and  the  complex  downwash  amplitude  is 

w  (x,y,z)  =  G0  exp  [  -i 2tt (z  +  za)/x  ]  (4) 

Introducing  the  reduced  frequency 

kr  =  ujC/2U  (5) 

and  eliminating  the  wave  length,  we  obtain 

w(x,y,z)  =  G0  exp  [  -i2kr(z  +  za)R/c  ]  (6) 

where 

R  =  U/Vg  (7) 
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We  next  consider  the  general  three-dimensional  problem.  If,  at  a 
particular  point  in  space,  the  gust  field  velocity  is  again  GU),  and 

—V 

the  velocity  vector  of  the  gust  front  is  Vg,  then  the  gust  vector  G  is 
G  =  (Vg/| Vg| )G(c)  (8) 

where  1 Vg 1  is  the  magnitude  of  Vg.  In  terms  of  the  direction  cosines 
of  the  gust  velocity  vector,  the  unit  gust  velocity  vector  is 

WN  (9a> 

=  icOSa  +  JCOSB  +  kcOSy  (9b) 

-*■ 

If  the  unit  normal  vector  to  a  lifting  surface  is  given  by  n,  then  the 
component  of  the  gust  vector  normal  to  the  surface,  i.e,  the  normalwash 
is 


Wg  =  n  •  G 

(10a) 

*  (n  •  Vg)G(c) 

(10b) 

■  e  G  (?) 

(10c) 

where  e  =  n  •  vg  is  a  parameter  that  depends  only  on  the  aircraft  geometry 
and  the  relative  angles  between  the  blast  and  the  aircraft. 


A  two-dimensional  variation  of  Figure  1  is  shown  in  Figure  2.  In 
the  three-dimensional  case,  the  penetration  distance  is  found  as  follows. 
Let  the  initial  vector  from  the  blast  origin  to  an  arbitrary  point  on  the 


aircraft  be 

r  =  T(x  +  xa)  +  j(y  +  ya)  +  k(z  +  za)  (11) 

and  the  velocity  vector  of  the  aircraft  be 

O'  =  -  Tu  (12) 

Then  the  penetration  distance  is  given  by 

5  =  Vgt  -  (Ut  +  r)  •  vg  (13a) 

=  (Vg  +  U  COSa)t  -  (x  +  Xa)cOSa 

-  (y  +  ya)cosB  -  (z  +  za)cosY  (13b) 


In  two-dimensions  6  =  90°  and  y  =  90°  -  a,  and  Equation  (13b)  leads  to 
the  result  shown  in  Figure  2. 

For  the  harmonic  case  we  again  assume  the  forms  in  Equation  (2). 

From  Equations  (2)  and  (13b),  we  find  the  frequency  to  be 

(0  =  2iT  (Vg  +  U  COSa)/A  (14) 

and  the  complex  normalwash  amplitude  to  be 

w  (x,y,z)  =  G0  e  exp  [  -i  2iri(x,y,z)A  ]  (15) 

where 

*(x,y,z)  =  (x  +  xa)  cosa  +  (y  +  ya)  cose  +  (z  +  za)  cosy  (16) 

Again,  we  may  introduce  the  reduced  frequency  kr  and  eliminate  the  wave 
length,  and  we  obtain 

w  (x,y,z)  =  G0  e  exp  (-  i  2  krtR/c)  (17) 

where 

R  =  U/ ( Vg  +  u  cosa)  (18) 
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3.  DOWNWASH  INDUCED  BY  VIBRATIONAL  MOTION 

The  normalwash  Induced  by  motion  is  found  from  a  substantial  deriva¬ 
tive  of  the  deflection  normal  to  the  surface 

w(x,y,z,t)  =  lif£  +  |£  (19) 

For  harmonic  motion 

h(x,y,z,t)  =  Re  [h(x,y ,z)exp{iut)]  (20) 

w(x,y,z,t)  =  Re  [w(x,y,z)exp(iat)]  (21) 

so  that  the  complex  amplitudes  of  the  normalwash  and  deflection  are 
related  by 

w(x,y,z)  =  U  ||  +  ijf  (22a) 

=  u(||*12krS)  (22b) 

Equation  (22b)  may  be  generalized  for  all  points  in  the  matrix  equation 


In  Equations  (19)  through  (23),  the  deflection  h  (and  h")  is  required  at 
the  collocation  points  appropriate  to  the  aerodynamic  theory,  in  the 
case  of  the  Doublet-Lattice  Method  at  the  three-quarter  chord  location  of 
each  aerodynamic  box,  and  in  the  case  of  Piston  Theory,  at  the  midchord 
of  each  box.  This  requirement  sets  a  further  requirement  for  interpolation 
from  the  structural  grid  points.  The  structural  deflections  are  given 
in  the  form  of  rigid  body  modes  and  vibration  modes  since  modal  methods 
are  usually  employed  to  reduce  the  number  of  degrees  of  freedom  in  a 
response  analysis.  The  deflections  of  the  structural  grid  points,  &,  are 
determined  from  a  series  of  rigid  body  and  vibration  modes  $n  as 

a(x,y,z)  =  rqn  *n(x,y,z)  (24) 

n 
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where  the  qn  are  the  amplitudes  of  the  vibration  modes.  In  matrix  form 
we  have 

W  »  M  jq}  (25) 
The  deflections  at  the  aerodynamic  collocation  points  are  found  by  inter¬ 
polating  among  the  structural  grid  points.  If  the  interpolation  function 
is  I(x.,y.,z.;  xk»yk»zk)  then  the  deflections  of  the  aerodynamic  points 


h(xj>yj>zj)  =  1  I(xj,yj,Zj;xk,yk,zk)A(xk,yk,zk)  (26a) 

k 

=  £  q^  E  I  ( x j ,yj  ,z j  ;x^ ,y^ ,z^)  4>n(xk»yk»zk)  (26b) 

or,  in  matrix  form 

|h}  =  [Ijk]M  (27a) 

=  Cljk3C^3  {q}  (27b) 


The  downwash  follows  by  substituting  Equation  (26b)  into  Equation  (22b) 
and  noting  that  the  streamwise  differentiation  only  involves  the  aerody¬ 
namic  coordinates  (the  j-set)  and  not  the  structural  coordinates  (the  k- 


set).  The  downwash  at  point  j  becomes 

Sj  a  u  l  %  l  (*f 1  i^) 


or,  in  general, 

'  |g|  .  u||ijk]  *  tijkl)  [♦]  111  (29) 

It  remains  to  discuss  the  form  of  a  suitable  Interpolation  function. 


A  linear  spline,  which  is  based  on  the  small  deflection  equation  of 
an  infinite  uniform  beam,  has  been  quite  useful  for  one-dimensional 
interpolation.  A  variety  of  two-dimensional  schemes  have  been  proposed 
over  many  years  all  of  which  have  had  one  shortcoming  or  another.  An 
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ingenious  resolution  of  this  classical  problem  in  Aeroelasticity  is  the 
surface  spline  of  Harder  and  Desmarais  (Ref.  9).  It  is  based  on  the  small 
deflection  equation  of  an  infinite  uniform  plate.  The  main  advantages 
of  the  surface  spline  are  that  the  coordinates  of  the  known  points  need 
not  be  located  in  a  rectangular  array  and  the  function  may  be  differentiated 
analytically  to  find  slopes. 

A  variation  of  the  linear  spline,  which  can  be  used  for  interpolation 
along  an  elastic  axis  is  developed  later.  It  utilizes  the  elastic  axis 
degrees  of  freedom;  i.e.,  the  normal  deflection,  twist  about  the  axis, 
and  the  bending  slope  of  the  axis,  in  order  to  interpolate  for  the  aero¬ 
dynamic  downwashes  in  the  streamwise  coordinate  system.  However,  we  will 
develop  the  surface  spline  equations  first. 

The  surface  spline  relates  the  deflections  Ak  at  the  structural  grid 
points  to  the  structural  loads  at  all  of  the  same  points. 

N 

Ak  =  a0  +  alXk  +  a2yk  +  r  Fa  Rkt  (30) 

£=1 

where 

Rki  =  V  *n  rk,  (31) 

and 

=  (xk  -  xj2  +  (yk  -  ya)2  (32) 

where  k  and  i  range  over  the  structural  deflection  points.  The  three 
additional  equations  to  determine  the  N  +  3  unknowns  are  the  stai'c 
equilibrium  equations. 


9.  Harder,  R.L.,  and  Desmarais,  R.N.  ."Intc'ipot'aWcn  Vitro  Spf<ra 

J.  Aircr.,  Vol .  9,  No.  2,  pp.  1 89-91 9,  19/2. 
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E  F^  =  0 

(33) 

Z  xfcF ^  =  0 

(34) 

1  =  0 

(35) 

In  matrix  form,  Equation  (30)  and  Equations  (33)  throuqh  (35)  may  be 
wri tten  as 

{A}=  [Pk]  {a}  +  [RkJ  |F|  (36) 

and 

[Pk]T  |F|  =  0  (37) 

where  [Pk]  denotes  the  planar  coefficient  matrix 

[Pk]  =  [  I  i  xk  •  yk]  (38) 

We  add  a  rigid  body  displacement  term  to  Equation  (37)  so  it  now  reads 

[B]  ja)  +  [Pk]T  jF}  =  0  (39) 

The  matrix  [B]  normally  should  be  zero,  except  in  the  case  where  yk  is 
constant,  as  for  a  fuselage,  an  external  store,  or  an  elastic  axis  beam; 
in  this  case  the  nonzero  value  of  [B]  chosen  as 


[B]  = 


0  0  0 

0  0  0 

.0  0  1 

avoids  a  singular  solution,  and  results  in  a  quasi-linear  spline  inter¬ 


(40) 


polation. 


Solving  Equations  (36)  and  (39)  simultaneously  yields 
M  =  -  ([B]  -  [Pk]T  [Rkd'1  [Pkl)'1  [Pk]T  [Ria]'1  |A}  (41a) 

-  da*]  lAl  (41b) 

{F|  =  [RkJ'1  (H-  Cpk]  ja})  (42a) 

-  [r kd"1  (  nj  -  [pk]  [iaj)  {a|  (42b) 

=  [I Fa]  w  (42c) 


*The  usual  linear-spline  is  a  series  of  continuous  cubics  and  differs  from 
the  degenerate  case  of  the  surface-spline.  Either  provides  an  acceptable 
interpolation  function. 
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The  deflections  at  the  aerodynamic  points  are  obtained  by  rewriting 
Equation  (36)  as 


{h}  =  [ P . ]  (a)  +  [R.S]  (F) 

J  J 

where  P.  are  the  planar  coefficients  for  the  aerodynamic  points 
J 

[P  ]  *  1  I  *  x.  *  y.  ] 

J  •  J  f  J 


and 


2  2 
d  =  r'  r 

jft  jit  jit 


(43) 

(44) 

(45) 


where 


rjt  ■  (xj  •  x«>2  *  (),j  • 1,1)2 


(46) 


where  "j"  ranges  over  the  aerodynamic  box  points  and  "it"  ranges  over  the 
structural  points. 

Substituting  Equations  (41b)  and  (42c)  into  Equation  (43)  and 
identifying  the  result  with  Equation  (27a)  leads  to  the  desired  inter¬ 
polation  matrix. 

"jk1  ■  "V  "as1  +  IV  "fa1  (47) 

The  derivative  of  Equation  (47)  is  also  needed  in  Equation  (29).  It  is 

(48) 


fef1 

il 

’  3P." 

3x1 

f-asl  1 

ffi 

f’Fsl 

L  J 

L  J 

L  J 

where 


and 


0] 


=  [2(x.  -  x  )(1  +  £n  r  )] 

J  *  j£ 

4.  AERODYNAMIC  INFLUENCE  COEFFICIENT  MODIFICATIONS 


(49) 


(50) 


An  analysis  of  the  blast  onset  pressure  field  is  presented  in  Appendix 
A.  The  gust  or  velocity  field  and  the  overpressure  field  are  different 
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aspects  of  the  blast  wave  and  are  not  independent.  To  be  completely 
accurate,  the  lift  induced  by  the  onset  pressure  field,  however  small 
it  is,  should  be  added  to  the  lift  induced  by  the  gust  field.  Presently, 
the  aerodynamics  module  calculates  only  the  loadings  induced  by  the  gust 
field.  When  it  becomes  desirable  to  include  the  additional  effects  of  the 
overpressure,  the  equations  of  Appendix  A  can  be  added  to  the  aerodynamics 
module. 

The  matter  of  lifting-  or  side-forces  on  bodies  with  elliptical  cross 
sections  has  also  been  the  subject  of  futher  investigation  since  the  pub¬ 
lication  of  Reference  4.  The  control  parameter  IBFS  determined  which  of 
two  methods  could  be  used  to  integrate  the  body  loads.  Experience  with  the 
option  IBFS=0  has  not  led  to  generally  satisfactory  results,  so  the  option 
I B FS = 1  for  circular  bodies  has  been  chosen  for  present  applications  with 
a  modification  for  elliptical  bodies  based  on  slender  body  theory.  The 
modification  for  elliptical  bodies  uses  the  circular  body  formula  for  lift, 
but  for  side-force,  multiplies  the  lift  by  the  square  of  the  cross  sectional 

aspect  ratio  (AR  =  height/width),  i.e.,  Fz  is  found  from  the  formula  in 

2 

Reference  4  for  IBFS  3  1,  and  F^  =  (AR)  F  The  effect  of  compressibility 
has  also  been  considered  in  a  modification  of  Slender  Body  Theory  using  the 
technique  of  Miles  (Ref.  10, p.  165).  Miles'  basic  equation  for  the  running 
load  3(F/q)/3x,  in  the  direction  of  the  upwash  w  (see  Figure  3),  is 

10.  Miles,  J.W.,  The.  Potential  Theory  of,  UnAt.ea.dij  Supe.uon.ic  Flow, 
Cambridge:  At  the  University  Press,  1959. 
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3(F/q)/3x  =  2^ 

(51) 

DM 

Dt  ‘ 

3M  .  jta)  y 

37  +  1U  M 

(52) 

M  =  1 

a  cos$  dt 

(53) 

c 


where  $  is  the  two-dimensional  potential  due  to  the  crossflow  w  and  is 
obtained  from  either  the  incompressible  or  compressible  crossflow  dif¬ 
ferential  equations.  Usually  the  incompressible  equation  is  used  for 
subsonic  flow;  however,  at  high  frequency  the  compressible  form  should  be 
used.  The  compressible  result  is: 


J  *  aw  S^(p)  cost 

(54) 

p  =  2kr  M  a/<f 

(55) 

Sl(p)  1  +  p$ 

(56) 

S  *  H^)(p)/H1(2)(p) 

(57a) 

=  <J0(p)  -  1Y0(p) >/(J1(p)  -  i (p) > 

(57b) 

in  which  "a"  Is  the  local  body  radius.  Placing  J  Into  the  expression 

for  M  gives 

M  =  a^w  S-j  (p )  tj> cos2t  dt 

(58a) 

*  w  a2w  S] (p) 

(58b) 

Then  substituting  M  Into  the  expression  for  the  running  load  gives 

2,r  Dt  ju  a  Sl<p)j 

(59a) 

■*iw 

(59b) 

in  which 

-  <>  -s2-  S/P) 

(60a) 

3p  _  2  kr 

3x  j-  dx 

(60b) 
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The  first  term  (with  S-j  =  1)  is  the  same  as  that  given  by  the  incompressible 
theory,  the  second  term  and  the  coefficient  S-j  are  the  compressibility 
corrections. 


First  order  piston  theory  for  motion  (Ref. 

pressure  coefficient  on  a  box  given  by 

.r  .  4  w 
ALP  ~  M  U 


8)  results  In  a  lifting 

(61) 


where  the  downwash  w  is  taken  as  the  value  at  the  midpoint  of  the  box 

centerline.  The  box  lift  Is  assumed  to  act  at  the  same  midpoint.  This 

assumption  has  sufficient  accuracy  because  at  high  frequency  the  downwash 

comes  primarily  from  the  plunging  motion  which  is  relatively  constant 

over  the  small  box  area.  However,  at  high  frequencies  there  is  substantial 

variation  of  the  gust  downwash  over  the  box  so  it  is  necessary  to  evaluate 

the  lift  by  integration  over  the  box  area.  The  assumotion  is  still  made 

that  the  lift  acts  at  the  box  center  and  is  justified  because  the  net  lift 

is  small  at  high  frequencies.  The  average  pressure  is  then  found  by 
dividing  the  integrated  lift,  by  tne  box  area.  The  geometry  of  a  typical 
box  is  shown  in  Figure  4. 

X1 
Yl 
x2 
c 


I  W 

Figure  3.  Cross  Flow  to  a  Body  Cross  Figure  4.  Lifting  Surface  Box  Geometry 
Secti on 
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aCp  = 


/'  4  w 

I  iijF 


dA 


(62) 


where 


Wg  =  Go0  exp  (-12  krn  R/c)  (63) 

and  the  Integration  Is  performed  over  the  box  area.  Carrying  out  the  Inte¬ 
gration  of  Equation  (62)  leads  to  the  following  expression  for  the  pressure 
coefficient. 


where 


4  G0e 

aCP  *  H  WT cosa  {(A31  C0St41  •  B31  sinTil 
-A42  cost£2  +  B42  sin^) 

+  1  (-B31  cost i]  -  A31  sinTi] 

+  B42  costs.2  +  A42  sinxZ2)} 

Ax  Is  the  average  box  chord  as  shown  in  Figure  4, 
t  =  2krR/c 

*m  =  xm  cosa  +  (y0  -  e  cosy)  cose  +  (z0  -  e  sinT)  cosY 

Amn  ■  (costAij^-l  )/TAJ.mn 
Bum  *  (slnTAfcmnJ/TAfcmn 


(64) 


and 


AS. 


•mn 


=  (*m  -  xn)  cosa  +  2e  (cosy  cose  +  sinY  cosy) 


(65) 

(66) 

(67) 

(68) 

(69) 


and  where  e  =  n  •  Vgj  i.e.,  the  component  of  the  gust  normal  to  the 
element. 


The  case  of  a  =  90°  (cosa  *  0.0)  requires  special  consideration. 
The  limiting  value  of  ACp  Is 

c„ .  e-<*0 1, „ 


VX1  -  X2  +  Xj  -  X, 


teg 


(70) 
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Substituting  these  expressions  into  Equation  (75)  gives 

t‘F/,)  ■  **«  5  C0Sa)  -  -  4)1  (so 

This  expression  gives  the  loading  at  a  particular  body  station.  Because 
of  the  variation  in  downwash  at  high  frequency  along  a  body  lifting  element, 
it  is  necessary  to  integrate  over  the  length  of  the  element  to  obtain  the 
total  load  to  be  used  in  the  generalized  force  calculation.  For  a  slender 
body  it  can  be  assumed  that  "a"  and  da/dx  are  constant  along  the  element 
and  the  average  values  are  taken  at  the  element  midchord  point.  Then 
the  integrated  force  on  the  element  of  length  ax  and  centered  at  x0  is 

ax/2 


AF/q  =  L  / 


■Ax/2 


Jdx 


where 


L  =  2air  li  (1‘Rc0So)  ■i2k7iM  2ff)l 


and 


x  =  x  -  x. 


The  integrated  force  becomes 

AF/q  _  o_.  w \  slnx 
~ AX  '  U(X  =  — 


where 


tt  (1-R  COSa) 

. /da 1_\| 

12kraMj)'  dx  2M  '\ 


(81) 

(82) 

(83) 

(84) 


T  =  Ax  krRcosa/c.  (85) 

When  "point"  running  loads  are  needed,  we  simply  set  (sint)/T  equal  to 
unity  in  this  formula. 


In  general,  an  arbitrary  burst  location  results  in  an  asymmetric 
loading  of  the  flight  vehicle.  Since  any  asymmetric  condition  can  be 
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found  as  a  linear  combination  of  a  symmetric  and  an  anti  symmetrical  con¬ 
dition  (assuming  linear  systems,  as  we  are),  both  the  symmetric  and 
antisymmetric  aerodynamic  analyses  need  to  be  performed.  In  Reference  4, 
either  is  optional  to  the  user,  but  for  the  present  applications,  the 
modifications  always  to  generate  both  have  been  added.  If  the  gust  normal - 
wash  on  the  right  side  is  wg  (x,y,z,7)  and  on  the  left  side  is  wg  (x,-y, 
z,-7),  then  the  symmetrical  normalwash  is 

Wgs  =  \  [wg(x,y,z,7)  +  wg(x,-y,z,-7)]  (86) 

and  the  anti  symmetrical  normalwash  is 

wga  =  \  Cwg(x,y,z,7)  -  wg(x,-y,z,-7)]  (87) 

5.  AERODYNAMIC  DATA  REQUIRED  FOR  VIBRA-6 

The  aerodynamic  data  required  from  Reference  4  for  the  blast  response 
analysis  consist  of  generalized  aerodynamic  forces  for  symmetric  and  anti¬ 
symmetric  modes  of  vibration  and  for  symmetric  and  antisymmetric  harmonic 
gust  fields  having  a  number  of  orientations.  The  orientations  depend  on 
the  burst  locations.  The  thirteen  burst  locations  are  shown  In  Figure  5 
and  the  standard  orientations  are  specified  by  the  direction  cosines 
tabulated  in  Table  1;  arbitrary  oritentations  must  be  specified  separately 
by  the  user.  The  generalized  forces  are  required  at  a  sufficient  number 
of  reduced  frequencies  that  Interpolation  will  lead  to  accurate  values 
at  Intermediate  frequencies. 

The  aerodynamic  geometry  is  also  required  In  order  to  perform  the 
integrations  that  lead  to  the  structural  responses.  That  is,  the  structural 
responses  are  required  in  the  form  of  shears,  moments  and  torques  at  selected 
critical  stations.  A  separate  stress  transformation  matrix  must  be  input  to 
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VI BRA-6 ,  for  the  determination  of  bendinp  and  shear  stresses  and  margins 
of  safety. 

6.  AERODYNAMIC  MONITORING  DATA 

The  computer  program  of  Reference  4  generates  aerodynamic  coefficients 
and  load  distributions  as  a  matter  of  course  for  each  mode  as  a  means  of 
checking  the  input  data.  This  is  necessary  since  generalized  forces  do 
not  have  enough  physical  significance  to  assist  the  user  in  detecting 
mistakes  in  input  data.  This  monitor  output  has  been  retained  In  the 
present  application  for  any  two  modes  selected  by  the  user;  e.g.,  sym¬ 
metric  pitching  and  antisymmetric  yawing.  The  monitor  output  for  the 
lifting  surfaces  includes  the  pressure  distributions,  spanwlse  distributions 
of  lift  coefficient,  moment  coefficient,  center  of  pressure,  and  loading 
(CjC/c),  and  total  lift  and  moment  coefficients.  The  monitor  output  for 
bodies  Includes  the  running  loads  (vertical  and  lateral)  and  the  total 
force  coefficients.  Finally,  the  total  force  and  moment  coefficients 
are  given  for  the  entire  aircraft.  The  specific  definitions  of  aerodynamic 
coefficients  are  presented  in  the  sample  output. 
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TABLE  1 

DIRECTION  COSINES  FOR  STANDARD  ORIENTATIONS 


ORIENTATION 

NUMBER 

Z* 

Zz 

1 

0.8660254 

0.0 

0.50 

2 

0.50 

0.0 

0.8660254 

3 

0 

0.0 

1.0 

4 

-0.7071068 

0.0 

0.7071068 

5 

-0.7071068 

0.0 

-0.7071068 

6 

0 

0.0 

-1  .0 

7 

0.50 

0.0 

-0.8660254 

8 

0.8660254 

0.0 

-0.50 

9 

0.7071068 

-0.7071068 

0.0 

10 

0 

-1 .0 

0.0 

11 

-0.7071068 

*  -0.7071068 

s 

0.0 

12 

0.0 

-0.7071068 

0.7071068 

13 

0.0 

-0.7071068 

-0.7071068 
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SECTION  III 
AIRCRAFT  IDEALIZATION 

Finite  element  methods  in  aerodynamics  (as  in  structures)  require  a 
compromise  between  accuracy  and  computing  time.  A  sufficient  number  of 
elements  must  be  used  to  obtain  acceptable  accuracy  but  not  so  many  as  to 
result  in  exorbitant  computational  costs.  Experience  has  been  gained 
with  the  Doubl et-Lattice/Image  Method  on  advanced  aircraft  configurations 
in  Reference  11  that  suggests  guidelines  to  achieve  this  compromise  in 
aoplications  to  gust  and  blast  response  analyses. 

The  basic  modelling  elements  include  lifting  surfaces,  slender  bodies, 
and  interference  surfaces.  These  are  consistent  with  the  assumptions  of 
linearized  aerodynamic  theory.  The  idealization  into  these  model  1 i ng 
elements  results  in  an  aircraft  geometry  that  is  much  simpler  than  the 
actual  geometry.  However,  the  essential  features  are  retained  such  as 
surface  planforms  and  body  volume  distribution,  although  some  secondary 
features  are  neglected  such  as  thickness  effects  of  the  lifting  surfaces. 
The  aerodynamic  coordinate  system  is  shown  in  Figure  6. 

I.  LIFTING  SURFACES 

Lifting  surfaces  are  idealized  as  plane  panels;  as  such  they  have  no 
thickness,  camber  or  twist,  but  do  have  dihedral.  The  effects  of  camber 
or  twist  and  angle  of  attack  are  accounted  for  in  the  linearized  smal1 
disturbance  boundary  condition  of  tangential  flow  at  the  surface;  i.e., 
the  so-called  normalwash  boundary  condition.  The  idealization  of  each 
lifting  surface  consists  of  small  trapezoidal  elements  (boxes)  arranged 
in  strips  parallel  to  the  freestream  (see  Figs.  7  and  8)  such  that 

II.  Giesing,  J.P.,  and  Kalman,  T.P.,  A eAodijnamic  Loadi  on  Advanced  Figkte’i 
AJjicfia^t,  Air  Force  Flight  Dynamics  Laboratory,  Report  No.  AFFDL-TR-73-45, 
October  1973. 

23 


1 


Figure  6.  Aerodynamic  Coordinate  System 


Figure  8.  Idealization  of  Lifting  Surface 


box  boundaries  lie  on  surface  edges,  fold  lines,  and  control  surface  edges. 

The  lift  on  each  box  is  represented  by  an  unsteady  horseshoe  vortex  having  its 
bound  leg  located  on  the  box  quarter-chord  and  its  trailers  extending  down¬ 
stream,  and  is  determined  by  satisfying  the  normalwash  boundary  condition  at 
a  control  point  at  the  three-quarter  chord  location  along  the  centerline  of 
each  box  (Fig.  7).  For  coplanar  or  near-coplanar  surfaces,  the  boxes  must  be 
aligned  in  the  streamwise  direction  so  that  the  control  points  are  centered 
between  the  trailing  vortices  (Fig.  8).  For  non-coplanar  surfaces,  the  strips 
need  not  be  aligned;  sufficient  perpendicular  separation  is  one  strip  width. 

For  intersecting  surfaces,  the  intersections  must  also  be  located  at  box 
edges  as  in  the  case  of  a  wing  pylon. 

The  lifting  surface  can  be  idealized  into  a  single  panel  if  there  are  no 
planform  (or  dihedral)  discontinuities;  otherwise,  the  planform  discontinui¬ 
ties  determine  the  subset  of  trapezoidal  panels.  A  panel  is  then  described  by 
the  coordinates  of  its  inboard  chord,  XI,  X2,  Yl,  and  Zl,  and  the  coordinates 
of  its  outboard  chord,  X3,  X4,  Y2,  and  Z2;  see  wing  panel  #1  in  Fig.  8.  The 
panel  then  has  a  dihedral Y«  arctan  [(Z2-Z1)/(Y2-Y1)  ]  ,  and  its  normal  vector 
is  n  =  -jsinY  +  tcosY,  where  3*  and  if  are  the  unit  vectors  in  the  y  and  z  direc¬ 
tions,  respectively.  The  panels  provide  the  means  of  generating  the  boxes 
within  the  computer;  the  boxes  are  generated  from  the  additional  input  informa¬ 
tion  of  the  spanwise  fractional  divisions  Tn  (n  =  1,  NS)  and  the  chordwise 
fractional  divisions  0i  (i  =  1,  NC).  The  panel  should  be  subdivided  so  that 
there  are  a  sufficient  number  of  boxes  per  wave  length;  i.e.,  the  box  chords 
should  be  approximately  less  than  c/6kr.  This  rule  was  relaxed  somewhat  for 
test  configuration  shown  on  page  53  to  reduce  cost.  Actually,  14  boxes  should 
be  used  on  the  root  strip  and  5  boxes  on  the  tip  strip  with  a  uniform  variation 
in  between.  This  corresponds  to  a  maximum  kr,  used  for  the  DLM,  of  2.0  based 
on  c  =  m.a.c.  The  minimum  number  of  boxes  on  a  chord  is  four.  The  number 
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of  strips  spanwise  should  be  chosen  so  the  box  aspect  ratios  are  of  order 
unity;  i.e.,  less  than  three  and  oreater  than  one-third. 

2.  SLENDER  CODIES 

The  idealization  of  fuselage,  nacelles,  or  external  stores  as  slender 
bodies  does  not  require  geometric  similarity  to  the  actual  body.  The 
slender  body  idealization  is  only  used  to  construct  the  boundary  conditions. 
It  is  recommended  that  the  best  elliptical  representation  of  the  body 
cross  sectional  area  be  used  independent  of  the  interference  body  geometry; 
this  idealization  is  permitted  to  pass  through  lifting  surfaces.  The 
elliptical  cross  section  is  specified  by  its  width  and  aspect  ratio  (ratio 
of  height  to  width).  The  cross  section  should  be  centered  on  a  straight 
line  along  the  planform  centerline  of  the  body;  the  centerline  coordinates 
are  at  YC  and  ZC  as  shown  in  Figure  9.  The  (small)  displacements  of  the 
body  during  motion  are  accounted  for  in  the  upwash  and  sidewash  boundary 
conditions. 

Slender  Body  Theory  is  used  to  calculate  the  axial  doublet  strength 
and  the  cross-flow  forces.  Hence,  the  maximum  loading  occurs  where  the 
cross-sectional  area  has  its  greatest  rate  of  change.  All  bodies, 
including  jet  engines,  must  be  idealized  as  having  pointed  noses.  This 
is  an  essential  assumption  of  Slender  Body  Theory.  However,  a  body  need 
not  be  closed  at  its  downstream  end  and  a  suitable  base  area  may  be  selected 
to  approximate  viscous  or  flow  separation  effects. 

In  addition  to  the  idealization  of  the  body  cross  sections,  the  body 
length  is  divided  into  segments.  The  segment  is  described  by  its  length 
and  its  cross-sectional  area  at  each  end.  The  cross-flow  force  is  assumed 
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to  act  at  the  center  of  the  segment  length.  The  segment  lengths  are  chosen 
from  considerations  of  streamwise  changes  in  geometry.  Smaller  lengths 
are  chosen  where  the  body  cross  section  is  changing  rapidly  and  in  regions 
of  maximum  interference  with  lifting  surfaces  as  near  the  wing  root  or  at 
a  nacelle-pylon  intersection. 

3.  INTERFERENCE  BODIES  AND  IMAGE  SYSTEMS 

A  cylindrical  tube  with  an  elliptical  cross  section  is  used  to  generate 
images.  All  attached  lifting  surfaces  are  connected  to  this  tube.  Within 
this  tube  is  placed  an  image  system  of  oscillatory  horseshoe  vortices. 

This  image  system  approximately  negates  the  flow  field  at  the  body  surface 
induced  by  lifting  surfaces.  The  images  of  lifting  surface  strips  which 
lie  close  to  the  interference  tube  are  of  prime  importance,  since  the  end¬ 
plating  effect  of  the  interference  tube  has  the  greatest  effect  on  these. 

However,  those  that  lie  further  than  two  interferance  tube  radii  can  usually 
be  neglected  in  the  interest  of  computational  efficiency.  However,  there 
should  be  at  least  two  adjacent  lifting  surface  strips  imaged  within  the 
tube.  For  this  reason  a  term  called  SCALER,  which  is  a  "non-dimensional 
image  radius"  is  an  input.  This  number,  when  multiplied  by  the  tube  radius, 
gives  a  radius  (from  the  tube  center  line)  within  which  lie  all  the  lifting 
surface  strips  which  possess  an  image  within  the  tube.  Strips  whose 
centerlines  lie  outside  of  this  radius  will  not  have  an  image  within  the 
tube  in  question. 

The  nondimensional  image  radius  is  not  input  per  body  but  is  the  same 
for  all  bodies.  Even  though  the  nondimensional  image  radius  (SCALER),  is 
the  same  for  all  bodies  the  actual  dimensional  image  radius,  a^  is 
different  for  each  body  since,  as  mentioned  above,  aj  is  the  result  of 


multiplying  the  term  SCALER  times  the  tube  radius  in  question.  Thus  the 
following  equation  holds: 

aj=  SCALER  a 

Thus  the  nondimensional  image  radius  (SCALER)  must  be  selected  such  that 
at  least  two  strips  of  the  adjacent  lifting  surfaces  be  imaged  for  each 
body/surface  combination.  If  this  is  not  a  problem  then  use  the  usual 
value  to  2.0. 

The  image  system  primarily  negates  the  effects  of  the  trailing  vortex 
system;  an  additional  slender  body  term  must  also  be  included  to  cancel 
the  residual  effects  arising  primarily  from  the  bound  vortex  system.  The 
interference  tube  is  segmented  into  elements  by  means  of  angular  divisions 
el  around  its  circumference  and  streamwise  divisions  along  its  length  as 
shown  in  Figure  10.  The  angular  positions  must  not  be  located  at  lifting 
surface/body  intersections,  and  there  must  be  at  least  two  positions;  four 
is  a  preferred  number  of  positions,  although  six  or  more  will  improve  the 
accuracy  of  the  average  upwash  and  sidewash  calculations  but  will  result 
in  excessive  computational  cost.  Since  the  interference  elements  increase 
the  computational  cost,  their  number  should  be  limited.  They  need  only  be 
used  near  lifting-surface/body  intersections  and  only  upstream  and  down¬ 
stream  within  a  short-length  of  the  lifting  surface.  In  situations  in  which 
high  accuracy  is  not  required;  eg.,  for  external  stores,  the  number  of 
interference  elements  might  be  severely  reduced. 

After  the  upwash  and  sidewash  have  been  averaged  around  the  inter¬ 
ference  body  circumference,  they  constitute  boundary  conditions  which  are 
applied  to  the  actual  slender  body  surface.  The  slender  body  cross- 
sectional  properties  must  be  input  (again)  at  the  interference 
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body  element  end  points.  This  input  is  called  the  interference  body 
radius  (or  width  for  elliptic  cross  sections)  distribution.  Once  the 
cross  flow  forces  on  the  slender-body/interference-body  combination 
are  calculated,  they  are  integrated  together  as  concentrated  forces 
at  the  slender  body  element  control  points  located  at  the  streamwise 
midpoint  of  each  segment. 

4.  ASSEMBLAGE  OF  BASIC  MODELLING  ELEMENTS  INTO  AN  AIRCRAFT  IDEALIZATION 

The  following  outline  summarizes  the  steps  necessary  to  idealize  a 
complete  aircraft  into  a  configuration  of  lifting  surfaces,  slender 
bodies,  and  interference  bodies. 

Step  1 .  The  aircraft  is  replaced  by  a  series  of  lifting  surface  panels, 
for  the  wings,  empennage  and  pylons,  and  interference  tubes  of  constant 
cross  section,  for  the  fuselage,  nacelles  and  stores.  An  example  is 
shown  in  Figure  11.  The  wing  is  made  up  of  three  panels.  The  first 
panel  is  inboard  of  the  pylon,  the  second  is  between  the  pylon  and 
the  aileron,  and  the  third  includes  the  aileron.  The  pylon  is  a  single 
panel  intersecting  the  wing  between  the  wing  first  and  second  panels. 
(Note:  it  is  not  necessary  that  a  pylon  intersection  be  between  wing 
panels,  but  only  at  a  strip  edge.)  The  horizontal  and  vertical  stabil¬ 
izers  are  each  a  single  panel  attached  to  the  fuselage  interference  tube. 

Some  liberty  must  be  taken  with  the  actual  configuration  in  order 
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to  accomplish  these  replacements.  For  example,  the  horizontal  tail  root 
is  moved  from  the  contour  that  fits  the  aft  fusel aqe  to  a  streamwise 
position  that  fits  the  fuselaqe  interference  tube.  This  may  actually  re¬ 
sult  in  moving  the  tail  outboard.  It  is  important  that  qaps  or  overlaps 
are  not  inadvertantly  permitted  between  adjacent  connected  components. 

The  nacelle  in  Figure  11  is  shown  pointed  at  the  nose  and  open  at  its 
base.  Bodies  should  be  idealized  with  pointed  noses  even  if  an  open  enqine 
nacelle  is  being  represented;  this  simulates  the  crosswise  flow  at  the 
intake,  bodies  need  not  be  closed  at  their  bases  and  an  open  base  can 
be  used  to  simulate  flow  separation. 

For  irregularly  shaped  bodies,  the  cross  sections  should  be  approxi¬ 
mated  by  ellipses  having  the  same  cross  sectional  area.  The  slender  body 
idealization  is  described  in  terms  of  the  body  width  and  aspect  ratio 
(the  constant  ratio  of  height  to  width  of  the  ellipse).  For  elliptic 
cross  sections  the  parameter  el  is  no  longer  the  anqle  at  which  the  pickup 
points  are  located.  It  is  a  parameter  such  that  the  pickup  point  locations 
are  given  as 

y  =  a  cosel 

z  =  a  AR  sinel 

Step  2.  Next,  the  surface  panels  are  divided  into  boxes.  Figure  11 
also  shows  these  subdivisions.  Four  to  six  chordwise  boxes  are  recommended 
but  if  a  panel  contains  a  control  surface  seven  to  twelve  are  advised. 

If  the  reduced  frequency  is  high,  a  box  chord  less  than  c/6  kr  is  necessary 
for  accuracy.  On  a  high  aspect  ratio  surface  10  to  15  spanwise  strips  are 
recommended.  For  low  aspect  ratio  tail  surfaces  at  least  8  strips  are 
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desirable.  One  or  two  strips  are  usually  sufficient  for  pylons  dependinq 
on  the  aspect  ratio.  In  general,  the  aspect  ratio  of  boxes  should  be  of 
order  unity;  e.g. ,  box  0.33  <  AR  <  3.0  is  suggested. 

At  surface  intersections,  e.g.,  on  a  T-tail,  the  leading  and  trailing 
edges  of  boxes  should  coincide.  However,  it  has  not  been  found  necessary 
to  align  box  end  points  on  adjacent  coplanar  panels,  although  it  probably 
is  desirable.  It  is  also  desirable  to  tailor  the  spanwise  strip  widths 
to  the  anticipated  spanwise  load  distribution;  narrower  strips  should  be 
used  in  regions  of  rapid  changes  in  loading,  e.g.,  near  the  tips.  Fiqure 
11  shows  the  horizontal  tail  divided  into  strips  aligned  with  the  wing 
strips.  This  is  a  requirement  for  a  coplanar  or  near-coplanar  tail  surface 
so  that  the  normalwash  collocation  points  on  the  tail  are  centered  between 
all  trailing  vortices.  If  vertical  separation  between  the  tail  surface  and 
the  wing  wake  is  greater  than  a  strip  width,  the  alignment  of  tail  strips 
with  those  of  the  wing  is  not  necessary. 

Step  3.  The  last  step  divides  the  bodies  into  slender-body  elements  and 
interference-body  elements.  Each  body  is  idealized  with  a  pointed  nose 
and  an  open  or  closed  base  and  with  an  elliptical  cross  section.  The 
elliptical  cross  section  has  a  variable  area  along  the  length  of  the  body 
which  is  equal  to  the  cross  sectional  area  of  the  actual  body  but  with  a 
constant  aspect  ratio  (height  to  width  ratio).  The  length  of  the  body 
is  divided  into  segments  with  shorter  segments  in  regions  of  maximum  rate 
of  change  of  cross  sectional  area  and  in  regions  of  a  liftinq-surface/body 
intersection.  A  sufficient  number  of  elements  are  required  to  describe 
accurately  the  body  cross  sectional  area  distribution.  However,  even  in 
regions  where  the  cross  section  is  not  changing,  a  number  of  elements 
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should  be  placed  to  establish  the  points  at  which  the  body  force  distribu¬ 
tion  is  calculated.  Figure  11  shows  a  number  of  elements  on  parts  of  the 
fuselage  and  parts  of  the  nacelle  where  the  cross  sections  are  constant. 

The  interference  tube  is  usually  divided  into  four  or  five  elements 
for  each  lifting-surface/body  Intersection.  Figure  11  shows  these  divisions 
by  a  dot-dash  line.  There  are  four  near  the  wing  and  four  near  the  tall. 
These  elements  compensate  for  the  residual  wing-body  interference  not 
accounted  for  by  the  image  system  in  cancelling  the  flow  through  the  fuse¬ 
lage  surface,  so  they  only  have  to  be  used  near  a  lifting-surface/body 
intersection.  In  cases  where  high  accuracy  is  not  required,  the  interference 
element  system  can  be  severely  reduced,  and  a  cost  savings  will  result.  The 
eioht  interference  elements  around  the  fuselage  (four  for  the  winq  and  four 
for  the  tail)  are  divided  around  the  circumference  at  three  points.  Only 
three  are  selected  so  that  the  computing  cost  is  reduced.  The  three  angular 
positions  around  the  tube,  el],  6l 2  and  8I3,  are  spaced  such  that  both  the 
horizontal  y-velocity  and  vertical  z-velocity  can  be  determined;  i.e., 
eli  and  91 2  9ive  the  horizontal  component  and  013  gives  the  vertical 
velocity.  If  there  were  no  vertical  fin  el^  =  90°  would  be  added  to  obtain 
a  more  accurate  average  vertical  velocity.  The  angular  positions  must  be 
chosen  away  from  lifting-surface  intersections  because  the  trailing  vortices 
at  these  locations  result  In  a  singular  calculation  of  crossflow  velocity. 

The  slender  body  radius  (or  width)  distribution  approximates  the  body 
whereas  the  Interference  tube  has  a  constant  cross  section.  The  slender 
body  radius  distribution  may  therefore  lie  within  the  tube  as  shown  in 
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Figure  11  or  It  may  pass  outside  the  tube  (not  shown),  and  It  may  even 
pass  through  lifting-surfaces.  A  total  of  27  slender  body  elements  on  the 
fuselage  and  10  on  the  nacelle  are  shown  in  Figure  11. 

5.  INTERPOLATION  OF  MODAL  DATA 

The  surface-spline  technique  described  In  Section  II,  Part  3,  is  used 
to  interpolate  for  the  aerodynamic  normalwash  distribution.  In  addition 
to  the  specification  of  aerodynamic  control  points,  which  has  been  discussed 
previously,  the  node  points  at  which  the  modal  deflections  are  given  (called 
aerodynamic  nodal  points)  must  also  be  specified.  Figure  12  shows  a  set 
of  nodal  points  for  the  configuration  of  Figure  11.  The  surface-spline 
fitting  through  the  deflections  at  the  nodes  to  obtain  deflections  and  slopes 
of  the  aerodynamic  boxes  is  illustrated  in  Figure  13.  It  should  be  empha¬ 
sized  that  the  surface-spline  is  very  accurate  for  interpolation,  but 
loses  accuracy  rapidly  as  extrapolation  distances  increase;  extrapolation 
should  be  avoided  whenever  possible. 

The  sign  convention  for  positive  deflections  on  bodies  is  upward 
(+z-d1 recti  on)  and  outboard  on  the  starboard  wing  (+y-d1rect1on).  On 
lifting  surfaces,  a  positive  deflection  is  normal  to  the  surface  where 
the  positive  normal  direction  is 

+n*  *  -jsiny  +  iTcosy 

in  which  the  positive  dihedral  y  is  given  in  Section  III,  Part  1. 

A  separate  surface-spline  is  fit  over  various  specified  areas;  the 
surface  fit  does  not  extend  beyond  the  boundaries  of  these  areas.  These 
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Superpanel  #2 


Deflection  of  Surface-Splines 


specified  areas  on  lifting  surfaces  are  referred  to  as  "superpanels"; 
on  bodies  they  are  called  "superbodies".  See  Figure  13.  In  Figure  12, 
Panels  No.  2  and  3  comprise  one  superpanel.  Dihedral  of  wing  panels  does 
affect  the  Interpolation  function  for  rigid  body  modes.  Thus,  panels  of 
dissimilar  dihedral  angles  are  not  combined  into  the  same  superpanel.  For 
panels  with  less  than  45°  dihedral,  the  surface-spline  is  based  on  the 
horizontal  (x,  y)  coordinates  of  the  vertically  projected  nodal  points; 
for  panels  with  greater  than  45°  dihedral,  the  spline  is  based  on  the  (x, 
z)  coordinates  of  the  horizontally  projected  nodal  points.  A  control 
surface.  Panel  No.  4  in  Figure  13,  must  always  be  a  separate  superpanel 
since  it  can  have  discontinuities  in  spanwise  deflections  and  streamwise 
slopes.  If  continuity  in  deflections  between  superpanel  edges  is  required, 
the  common  nodal  points  must  be  specified  twice,  once  on  each  superpanel. 

An  example  is  the  control  surface  hinged  at  its  leading  edge  (see  Fig. 

13);  the  nodal  points  on  the  hinge-line  must  be  specified  once  on  the  wing 
superpanel  and  again  on  the  control  surface  superpanel.  A  second  example 
is  the  common  side  edge  between  superpanels  #1  and  #2. 

A  second  type  of  interpolation  procedure  for  superpanels  is  obtained 
using  the  elastic  axis  degrees  of  freedom,  h,  the  normal  bending  deflection 
of  the  elastic  axis,  a,  the  twist  about  the  axis,  and  e,  the  bending 
slope  of  the  axis.  These  degrees  of  freedom  are  illustrated  in  Figure  14. 
These  elastic  deflections  are  given  at  nodal  points  along  the  elastic  axis, 
but  not  at  the  end  points  (XA,  YA,  ZA)  and  (XB,  YB,  ZB). 

As  with  the  spline  Interpolation,  control  surfaces  must  be  defined  as 
separate  superpanels.  If  absolute  (rather  than  relative)  control  surface 
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deflections  are  used,  the  wing  and  control  surface  superpanels  are  mutually 
exclusive  as  before.  However,  if  relative  control  surface  deflections  are 
used,  the  control  surface  superpanel  is  superimposed  on  the  wing  superpanel. 

In  this  case,  the  aerodynamic  boxes  on  the  control  surface  lie  on  both 
superpanels.  In  the  input  procedure  a  selection  is  made  among  the  degrees 
of  freedom  (h,  a,  and/or  e)  that  are  to  be  used;  e.g.,  if  relative  coordi¬ 
nates  are  used  for  a  control  surface,  then  only  a  might  be  used.  The  nodal 
points  along  the  elastic  axis  are  normally  the  mass  reference  points. 

However,  for  control  surfaces,  the  computer  program  automatically  calculates 
nodal  points  when  the  superpanel  is  flagged  as  "associated"  with  another 
superpanel  since  there  are  no  mass  reference  points  on  the  control  surfaces. 

Fiaure  14  gives  a  graphical  illustration  of  how  the  control  surface  nodes 
are  determined.  Lines  are  drawn  normally  to  the  elastic  axis  of  the  associ¬ 
ated  superpanel  and  not  along  radial  lines.  The  points  at  which  they  intersect 
the  control  surface  superpanel  elastic  axis  are  the  generated  nodes. 

Each  superpanel,  whether  it  be  a  wing,  control  surface,  tail,  etc., 
is  divided  into  "sections".  A  "section"  Is  an  area  of  a  superpanel  in 
which  a  specific  way  of  generating  the  wing  surface  deflection  and  slope 
(in  the  x-direction)  Is  employed.  Usually,  the  deflection  of  the  super¬ 
panel  is  given  by  line  generators  normal  to  the  elastic  axis  (see  Figure 
14,  superpanel  section  #2).  However,  special  provision  has  been  made  for 
the  root  and  tip  areas.  In  these  areas  a  radial  or  fan  type  of  line 
generator  is  used  so  that  the  line  generator  lies  parallel  to  the  Inboard 
and  outboard  chords  of  the  "section"  with  a  monotonic  variation  In  between. 
Figure  14  Illustrates  the  radial  line  generators  in  superpanel  sections  #1 
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and  #3.  In  this  case  the  outboard  chord  of  section  #1  is  normal  to  the 
elastic  axis  as  is  the  inboard  chord  of  section  #3.  However,  this  Is  not 
necessary. 


In  order  to  find  the  deflection  and  slope  at  all  points  on  the  surface, 
an  interpolation  procedure  must  be  applied  to  h,  a,  and  e  along  the  elastic 
axis.  A  modified  one-dimensional  spline  is  used  for  this  purpose.  This 
type  of  curve-fitting  function  is  not  adequate  for  extrapolation.  To 
overcome  this  difficulty,  a  linear  extrapolation  from  the  first  two  and 
last  two  nodes  is  made  to  obtain  the  deflection  at  the  elastic  axis  end 
points.  If  { h }  is  the  set  of  N  bending  deflections  at  the  nodes,  excluding 
the  elastic  axis  end  points,  then  the  set  of  (N  +  2)  nodal  deflections  (h"} 
which  include  the  end  points  may  be  written 


{h}  =  [T]  {h} 


where  the  transformation  matrix  [T]  appear  as 


a  b  0 


[T]- 


1  0  0 
0  1  0 
0  0  1 


0  0 
0  0 
0  0 
0  0 


(88) 


(89) 


i 


i 


0  0  0.  .  .  0  1 
0  0  0.  .  .  c  d 
in  which 


a  * 

1  -  b 

(90) 

c  * 

1  -  d 

(91) 

b  * 

(T1  -  ^/(tj  "  t2^ 

(92) 

c  ■ 

(tn  +  2  -  tn)/<tn  +  i  ‘  TrP 

(93) 
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where  the  -r-j  are  the  nodal  coordinates,  including  the  end  points  as 
measured  along  the  elastic  axis  from  XA,  YA,  ZA  (see  Figure  14);  i.e., 
f|  =  0.  The  same  transformation  applies  to  obtain  the  slopes,  (a)  from 
{a}, and  {?}  from  { e } ,  so  we  write 

{a}  =  [T]  {a}  (94) 

{i}  =  [T]  (e)  (95) 

Separate  one-dimensional  splines  are  then  fit  through  the  values  of 
h,  a,  and  e  at  all  of  the  node  points  and  the  end  points  of  the  elastic 
axis,  as  illustrated  in  Figure  14.  From  the  variation  of  h,  a,  and  e 
along  the  elastic  axis,  which  we  denote  as  h(i),  a(t),  and  e(-r),  we  may 
find  the  deflection  H  of  any  point  on  the  superpanel  surface  and  its 
streamwise  slope  dH/dx.  There  are  two  cases  depending  on  whether  the  point 
is  on  a  parallel  line  generator  normal  to  the  elastic  axis  or  on  a  radial 
line  generator.  For  points  on  the  parallel  line  generators,  the  deflec¬ 
tion  is 

H  =  h(-r)  +  y  sinA  -  (x-XA)  cosa  a(x)  (96) 

where  x  and  y  are  given  in  Equations  103  and  104. 

The  streamwise  slope  is 

=  -e(i)s1nA  -  a( t )cosa  +  sinA  (y  sinA  -  (x  -  xA)  cosa)  (97) 
where  A  is  the  sweep  angle  of  the  elastic  axis  In  the  plane  of  the  surface 
and  t  is  the  distance  along  the  axis  from  the  origin  at  (XA,  YA,  ZA). 

For  points  with  radial  line  generators  from  an  origin  at  (X0,  Y0,  Z0), 
the  deflection  Is 

H  *  h(x)  -  (p  -  Pb)cos  (A  +  4)  [a (t )  +  e(i)  tan(A  +  $)]  (98) 
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and  the  streamwise  slope  is 


dH 

dx 


where 


ob  = 
and 

x  = 

y  = 
xo  = 

YO  = 

Po  " 

*0  = 
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y  is  the  dihedral  angle  of  the  plane  containing  the  elastic 
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SECTION  IV 


APPLICATIONS  TO  RIGID  AIPCRAFT 

A  number  of  applications  have  been  made  to  a  two-dimensional  airfoil  and 
to  a  typical  large  twin-enqine  transport  that  is  assumed  to  be  rigid.  The 
applications  are  intended  to  exhibit  the  convergence  characteristics  of  the 
Fourier  transformation  and  the  effects  of  the  transition  from  the 
Doublet-Lattice  Method  (DLM)  to  Piston  Theory  (PT). 

The  number  of  calculated  aerodynamic  solutions  depends  on  the  blast 
orientation.  An  estimate  of  the  required  reduced  frequency  spacino  is 

Akr  =  irU+M^cosaJ/fcjMj/c  )  (111) 

where  l  is  the  slant  distance  from  the  pressure  point  in  Question  to  the 
blast  plane  as  it  passes  through  the  aircraft  origin.  To  render  all  pressure 
solutions  accurate,  the  maximum  value  of  J  is  usually  used  in  this  formula. 
Also,  Afj  is  the  number  of  solutions  per  cycle,  recommended  in  the  range  of 
6  to  8.  The  maximum  frequency  that  must  be  considered  is 


/l  +  M  cosa\  /l  +  McoscA 

\  M  c/U  )  orkr=  M  ) 

\  »  °°  /  max  \  co  / 


(112) 


where  ^  is  recommended  as  2  or  3. 

The  transition  gap  between  the  two  aerodynamic  theories  results  in 
differing  frequency  response  curves.  The  smaller  the  gap,  the  more  pronounced 
is  the  difference  between  the  frequency  responses  and  the  larger  are  the 
ripples  in  the  time  domain. 

Calculated  results  are  shown  in  Figures  15  through  31.  Each  solution 
is  based  on  a  step  function  escitation  frequency  response  because  the  step 
function  is  a  good  approximation  to  the  frequency  response  of  a  blast  wave 
(see  Figure  19)  and  the  time  histories  will  approach  steadystate  values 
that  are  known. 
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Figures  15  through  17  show  results  for  a  two-dimensional  airfoil. 

Figure  15  shows  the  convergence  of  the  DLM  to  the  PT  for  pitching,  plunging, 
and  gusts  at  45°  and  90°  angle  of  attack  in  the  frequency  domain. 

The  term  f  ,  in  many  of  the  figures  to  follow,  is  the  frequency,  in 
max 

cps  or  Hz,  at  which  the  frequency  response  is  terminated.  Such  a  termina¬ 
tion  leads  to  a  small  error  originating  in  the  inverse  Fourier  Transforma¬ 
tion.  This  error  becomes  smaller  as  f  is  increased.  Figure  16(a)  corn- 

max 

pares  the  transient  lift  (with  effectively  infinite  fmax  in  the  PT)  to  re¬ 
sults  obtained  by  Lomax,  etal  .  (Ref.  12).  The  calculations  followed 
Zartarian's  suggestion  (Ref.  7)  of  obtaining  the  PT  contribution  analyti¬ 
cally  and  the  contribution  from  the  difference  between  the  DLM  and  PT  from 
the  numerical  Fourier  transform.  Figure  16(b)  shows  the  pressure  at  two 
points  near  the  leading  and  trailing  edges  induced  by  a  45°  gust.  The 
irregularities  in  the  curves  show  the  effects  of  the  DLM-PT  transition  at 
=  5.0.  The  agreement  with  PT  is  seen  at  the  appropriate  intercept  times 
for  each  point.  Figure  17  presents  a  convergency  study  as  a  function  of 
maximum  kr  with  a  vertical  gust.  The  DLM-PT  transition  is  at  kr  =  10.0. 

Slow  convergence  is  seen  in  this  case  at  high  subsonic  Mach  number  (M  =  0.8). 

Figure  18  gives  the  idealization  of  a  three-dimensional  transport;  the 
test  case  in  Volume  I,  which  is  being  used  to  gain  insight  and  experience  in 
the  use  of  the  aerodynamic  aspects  of  VIBRA-6.  Even  though  the  idealization 
is  crude,  much  valuable  information  has  been  gained.  The  dots  shown  on  the 
transport  idealization  are  the  points  at  which  pressures  or  running  loads 
are  calculated  on  the  rigid  vehicle.  No  problem  has  been  found  on  the 
elastic  vehicle  once  the  aerodynamic  frequency  response  has  been  well 


12.  Lomax,  H.  ,  Heaslet,  M.A.,  and  Sluder,  L.,  "Tkc  Indiciat  Li  fit  and  Pitch¬ 
ing  Moment  fioK  a  Sinking  ok  Pitchinq  Two-Vimc.ni  ionai  Wing  PCijing  at  Subionic 
ok  SupeKionic  Spcedi",  NACA  TN  2403',  July  1951. 
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defined  for  the  rigid  vehicle.  Notice  the  row  of  dots  on  the  section 
located  at  Y/(b/2)  =  0.71.  This  allows  for  the  integration  of  a  section 
c*  at  this  spanwise  station. 

Figure  19  provides  the  basis  for  using  the  step  function  blast  approxi¬ 
mation.  A  typical  blast  wave  does  not  differ  significantly  from  the  step 
function  (notice  that  it  approaches  the  step  function  very  quickly;  i.e., 
before  f  =  5.0  [kr  =0.42])  and  the  step  function  leads  to  known  asymptotic 
steady-state  values.  The  term  g  is  the  modulus  of  the  frequency  response 
of  a  typical  neuclear  blast  gust  velocity. 

Figure  20  shows  the  frequency  response  of  a  leading  edge  pressure  from 
a  vertical  blast.  The  frequency  response  curve  is  obtained  by  fitting  a 
spline  curve  through  specific  hard  points  (in  frequency).  The  frequency 
hard  points  are  shown  as  dots  along  the  abscissa.  The  transition  of  this 
curve  between  the  Doublet  Lattice  and  Piston  Theories  is  simply  obtained 
by  leaving  a  gap  in  hard  point  frequencies  between  the  two  and  allowing  the 
spline  to  fill  this  gap.  Different  curves  will  result  for  different  gaps. 
Notice  in  Figure  20  the  two  different  curves  in  the  gap  regions  for  the 
real  part.  Two  separate  curves  are  not  noticeable  for  the  imaginary  part. 

The  effects  of  these  gaps  are  seen  in  subsequent  figures.  In  general, 
the  smaller  the  gap  is,  the  greater  is  the  discontinuity  in  the  frequency 
response,  and  the  irregularities  in  the  time  history  are  more  extensive. 

The  frequency  response  of  a  leading  edge  pressure  is  shown  in  Figure  21. 
At  low  frequency  (f  <  10  Hz),  the  original  distribution  of  calculated  points 
did  not  contain  a  sufficient  number  (see  dots)  to  be  able  to  predict  the 
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peaks  accurately  by  interpolation.  Subsequently,  three  new  DLM  points  were 
added  signified  by  x's  on  the  axis.  In  the  PT  region,  the  same  problem  is 
seen  in  that  the  interpolated  curve  is  not  a  good  sinusoid.  No  new  points 
were  added  however  since  convergence  of  the  inverse  Fourier  transformation 
occurred  at  fmax  =  24  Hz  and  for  this  blast  orientation  piston  theory  is  not 
needed.  The  effect  of  the  insufficiency  of  calculated  points  is  seen  in 
later  figures. 

Figure  22  shows  excellent  results  for  a  60°  blast.  There  is  a  suffi¬ 
cient  number  of  calculated  points  and  the  DLM-PT  transition  has  no  signifi¬ 
cant  discontinuities. 

Results  in  Figure  23  show  the  lift  time  history  for  various  fmax 
values;  i.e.,  fmax  =  24,  120  and  300  Hz.  All  curves  lie  on  top  of  one 
another  which  shows  the  case  of  fmax  -  24  Hz  was  converged.  The  precursor 
in  the  lift  before  blast  interception  at  the  section  deserves  further 
investigation.  It  may  be  a  property  of  the  numerical  Fourier  transform. 

In  Figure  24,  we  see  another  converged  result  for  fmax  =  120  for  the 
case  of  a  60°  head  on  step-blast  gust.  The  transition  frequency  is  at 
f  =  24  Hz  as  in  Figure  22. 

Figure  25  presents  another  convergence  study  for  a  vertical  blast. 

The  figure  shows  that  fmax  is  an  important  parameter  for  the  vertical  blast 
In  determining  the  initial  pressure  rise.  This  is  a  general  property  of 
Fourier  transforms;  i.e.,  the  high  frequency  behavior  determines  the  initial 
response,  and  the  low  frequency  behavior  determines  the  approach  to  the 
steady-state  response.  This  figure  shows  that  the  three-dimensional  case 
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converges  faster  in  fmax  than  does  the  two-dimensional  case.  This  case 
has  converged  at  f_.„  =  300  (k_  =  25)  whereas  the  two-dimensional  case  of 

max  r 

Figure  17  converges  somewhere  between  kr  =  25  and  75. 

Figures  26  through  31  show  plots  obtained  with  the  computer  printout. 
Figures  26  through  23  present  pressures  at  three  different  points  on  the 
wing.  Figures  29  and  30  give  loadings  at  two  fuselage  stations,  and 
Figure  31  shows  pressures  on  the  horizontal  tail.  Figure  26  shows  pressure 
time  histories  for  an  outboard  point  forward  of  the  midchord  (x/c  =  0.406, 
2y/b  =  0.71)  for  three  blast  orientations.  The  vertical  line  shows  the 
•'ntercept  time  for  each  blast  to  reach  the  point.  The  initial  pressure 
rise  can  be  compared  to  the  Piston  Theory  value  shown  and  the  approach  to 
the  steady-state  solution  is  seen  at  large  values  of  time.  The  pressures 
are  negative  in  Figure  26{a)  because  the  blast  is  from  overhead. 

Figure  27  shows  pressures  at  the  same  spanwise  station  but  farther  aft 
(x/c  =  0.657,  2y/b  =  0.71)  for  two  blast  orientations.  Reasonable  agree¬ 
ment  with  both  Piston  Theory  and  the  steady  state  solution  is  observed. 
Figure  27(b)  also  shows  the  effect  of  increasing  the  number  of  DLM  points 
in  the  frequency  response  plot  (Figure  21).  The  three  additional  DLM 
solutions  improve  the  time-history  although  more  improvement  is  still 
needed. 

Figure  28  shows  pressures  at  the  same  outboard  station  but  at  a  point  near 
the  trailing  edge  (x/c  =  0.906,  2y/b  =  0.71)  for  an  overhead  blast.  The 
figure  shows  the  effect  of  two  transitions  (gaps)  between  DLM  and  Piston 
Theory.  The  transition  begins  at  24  Hz  in  both  cases  but  a  gradual  transi¬ 
tion  (large  gap)  reduces  the  irregularities  over  the  abrupt  transition 
(small  gap). 
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Figure  29  shows  the  running  load  at  the  fuselage  nose  for  three  blast 
orientations.  Figure  29 (a )  shows  a  large  effect  of  gap  size.  A  moderate 
gap  reduces  the  variations  in  the  time  history. 

Figure  30  gives  the  fuselage  loading  toward  the  aft  end  for  one  blast 
orientation.  A  similar  effect  of  gap  size  is  seen  again. 

Figure  31  shows  the  pressures  at  a  point  of  the  horizontal  tail  near  the 
leading  edge  and  tip  (x/c  =  0.0625,  2y/b  =  0.888)  for  three  blast  orienta¬ 
tions.  The  curves  in  Figures  31(a)  and  (b)  are  seen  to  be  well  behaved. 

The  pressures  begin  to  increase  to  the  steady  state  values  for  the  isolated 
tail  until  the  wing  downwash  reaches  the  tail  at  which  time  the  pressures 
decrease  to  the  steady  state  wing-tail  interference  value  of  approximately 
one-half  the  value  for  the  isolated  tail.  The  effect  of  adding  three  DLM 
points  is  shown  in  the  solid  curve  in  Figure  31(c).  The  general  effect 
is  to  smooth  the  curve  but  the  remaining  irregularities  are  unexplained. 
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DLM  (2D)  Theory 
Piston  Theory 


Figure  15.  Frequency  Response  of  a  Two- Di mens iona 
Airfoil  in  Pitch  and  Plunge,  and  in 
Two  Blast  Gust  Fields 


Figure  16  Time  Histories  of  (a)  Lift  and  (b)  Pressures,  Due  to  Airfoils  Encountering  Sharp  Edgti 
Gusts,  (a)  =  0.8,  a.  =  90  Deg  and  (b)  =  0.5,  o„  -  45  Deg. 


_  NACA  TN2403 

-  Present  Method 

jkr  <  10  DLM  Theory 
)kr  >  10  Piston  Theory 
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Figure  17. 


Time  History  of  Lift  Due  to  a  Step- 
Function  Blast  Gust  at  90°. 


•  Time  Histories 

•  Calculated  at  Dots 


Figure  18.  Doublet  Lattice  Method  (DLM)  Aerodynamic 
Idealization.  Also  Points  at  Which  Time 
Histories  are  Given,  a  re  Illustrated. 
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Figure  19.  Frequence  Response  of  a  Nuclear 
Blast  Gust  Velocity  Normalized 
by  that  Due  to  a  Step-Function. 
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Figure  23.  Section  Lift  Coefficient  Time  History  Due  to  a  45°  Tail-On  Step  Blast  Gust 


300,  120  Hz 


Section  lift  Coefficient  Time  History  Due  to  a  60°  Nose-On  Step  Blast  Gust 


25.  Section  Lift  Coefficient  Time  History 
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Figure  26.  Pressure  Time  Histories  for  a  Point  Located  on  the 
Wing  at  x/c  «  0.406  and  Y/(b/2)  *  0.71  (Station  7) 
for  Various  Blast  Orientations.  M»  ■  0.8. 
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Figure  27.  Pressure  Time  Histories  for  a  Point  Located  on  the 
Wing  at  x/c  -  0.657  and  Y/(b/2)  -  0.71  (Station  9) 
for  Various  Blast  Orientations.  M«  ■  0.8 
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Figure  29.  Time  History  of  the  Running  Load  at  the  Fuselage  Nose 
(Station  12)  for  Various  Blast  Orientations.  ■  0.8 
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Figure  31.  Pressure  Time  Histories  for  a  Point  Located  on  the 
Horizontal  Tall  (Station  3)  for  Various  Blast 
Orientations.  M..  «  0.8 
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SECTION  V 
PROGRAM  INPUT 

1.  DESCRIPTION  OF  INPUT  DATA 

The  program  input  consists  of  the  following  five  groups  of  data: 
general  data,  panel  data,  body  data,  modal  spline  interpolation  data,  and 
gust  data.  The  general  data,  panel  data,  and  modal  spline  interpolation 
data  must  always  be  input,  but  the  body  data  and  gust  data  are  optional 
as  specified  by  flags  in  the  general  data. 

The  general  data  consist  of  flags  for  dimensioning  and  program 
control,  locations  of  and  deflections  at  the  aerodynamic  nodal  points, 
and  constants.  The  panel  and  body  data  consist  of  the  input  necessary  to 
define  completely  each  panel  and  body.  The  modal  spline  interpolation 
data  consist  of  the  data  to  relate  the  aerodynamic  nodal  points  to  the 
panel  and  body  data.  The  gust  data  allow  nonstandard  blast  orientations 
to  be  input  for  the  gust  calculations. 

A  detailed  description  of  the  input  data  is  given  in  Volume  I  of 
this  report. 

2.  TEST  CASE  INPUT  DATA 

The  configuration  discussed  in  Section  III  (Figure  11)  is  used  as 
the  basis  for  the  test  case  presented  here.  In  the  interest  of  simplifying 
the  input  and  saving  space  for  the  output  this  configuration  has  been 
reidealized  with  fewer  lifting  surface  boxes,  fewer  body  elements  and  the 
tail  surfaces  missing.  Also  fewer  nodal  points  (at  which  modal  deflections 
are  given)  are  used.  Figure  32  presents  an  illustration  of  this  ideali¬ 
zation.  Table  2  tabulates  the  important  quantities  required  for  input 
while  Table  3  presents  the  actual  input  data. 
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TABLE  2 


TEST  CASE  DATA 


General  Data 

Mach  No.  =  0.85 
Ref.  Chord  c  =  2.0 
Ref.  Area  A  =  9.0 
Reduced  Frequency  kr  =  0.1 
Ref.  Semi  span  S  =  3.8 

Wing  Data 
Panel  #1 

* 

leading  edge  inboard  XI  =  0,  Y1  =  .8,  Z1  =  .6 
trailing  edge  inboard  X2  =  2.0 
leading  edge  outboard  X3  =  .6,  Y2  =  2.3,  Z2  =  .6 
trailing  edge  outboard  X4  =  2.1 

Chordwise  divisions,  in  fractions  of  chord,  e  =  0,  0.5,  1.0 
Spanwise  divisions  in  fractions  of  span,  t  =  0,  0.45,  1.0 

Panel  n 

XI  =  .6,  X2  =  2.1 ,  X3  =  1.2,  X4  =  2.2 
Y1  =  2.3,  Y2  =  3.8,  Z1  =  .6,  12  =  1.2 
9  =  0,  0.5  1.0;  t  =  0,  0.5,  0.75,  1.0 


Panel  #3 

XI  =  .6,  X2  =  2.1,  X3  =  0,  X4  =  1.5 
Y1  =  2.3,  Y2  =  2.3,  Z1  =  .6,  Z2  =  0 
9  =  0,  0.5  1 .0;  t  =  0,  1 .0 


*  Specific  units  not  required; however  they  all  must  be  consistent. 
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Body  Data 


Fuselage 

fuselage  center  coordinates  Y  =  0,  1  =  0,  tube  radius  a  =  1.0 

C*  V 

slender  body  element  endpoints  XIS  and  radii  RS  and  interference  element 
endpoints  XII  and  radii  RI 


XIS 

RS 

XII 

RI 

-2 

0 

-1 

1.0 

-1.8 

.65 

0 

1.0 

-1.4 

.95 

1 

1.0 

-1.0 

1.0 

2 

1.0 

0 

1.0 

3 

1.0 

1.0 

1.0 

2.0 

1.0 

3.0 

1.0 

4.0 

0.5 

5.0 

0 

Angular 

flow 

field  points 

el  = 

Aspect  ratio 

of  cross-section 

Nacelle 

Yc 

"  2'3’  Zc  = 

-.5, 

XIS 

RS 

XII 

RI 

-.4 

0 

0 

.5 

-.2 

.4 

.75 

.5 

0 

.5 

1.5 

.5 

.75 

.5 

1.5 

.5 

01  =  45°, 

,  135C 

\  225°,  315° 

Aspect  ratio  of  cross-section  =  1.0 
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TABLE  3 


TEST  CASE  INPUT  DATA 


JPOUF_NOi_l  General  Data 

ITEM_N01__.l  Data  Designator  Card 

Aiio 

ITEM  NOS.  2-4  Dimensioning  Data 

14  2  1  1  3 

3  6  3  2  12 

2  13  6 

ITEM  NO.  5  Aerodynamic  nodal  point  coordinates  x,y,z 


0.  0 

0.  8 

0.6 

2.  0 

0.8 

0.  6 

0.  6 

2.  3 

0.  6 

2.  1 

2.3 

0.6 

0.  6 

2.3 

0.6 

2.  1 

2.3 

0.  6 

1.  2 

3.8 

1.2 

2.  2 

3.8 

1.  2 

0.  3 

2.  3 

0.3 

1.8 

2.3 

0.  3 

2.  0 

0.0 

0.0 

5.  0 

0.0 

0.0 

0.  0 

2.  3 

-0.5 

1.5 

2.  3 

-0.5 

I1£M_N01_6  Modal  Deflection  data  PHINA,  PHIZA,  PHI YA 

Mode  1 

0.  0 

-2.  0 
-0,  6 
-2.  1 
-0.5571 
-1.  950 
-1.114 
-2.  043 
0.  0 
0.  0 

2.0  0.0 

-5.0  0.0 

0.0  0.0 

-1.5  0.0 

Mode  2 

0.  0 
0.  0 
0.  0 
0.  0 
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0.  0 
0.  0 
0.  0 
0.  0 
0.  3 
1.  8 


0.0 

0.0 

0.0 

0.  0 

0.0 

0.  0 

0.0 

1.  5 

Mode  3 

0.  0 
0.  0 
0.  0 
0.  0 

-0. 2228 
-0. 7799 
-0. 4456 
-0.  8170 

0.  3 

1.8 


0.0 

-2.  0 

0.  0 

5.0 

0.0 

0.0 

0.0 

1.  5 

ITEM  NOS,  7-8  Control  Flags 


0 

0 

0 

0 

1 

0 

0 

0 

ITEM  NO.  9  Case 

Data,  Mach 

number.  A, 

S,C, XM, SCALER 

0.  85  9.0 

3.8 

2.0 

0.0 

ITEM  NO.  10  Reduced  Frequencies 

0.  1 

2£Oup_no1__2 

Panel 

Data 

ITEM  NOS.  1-5  XI 

•  •  • Z  2  §  NC  § 

NS,  IGPUP, 

Theta,  Tau 

Panel 

1 

0.0  2.0 

0.6 

2.  1 

0.8  2.3 

0.6 

0.6 

2 

2 

1 

0.0  0.5 

1.0 

0.0  0.45 

1.0 

Panel 

2 

.6  2.1 

1.  2 

2.  2 

..  }  3.8 

0.6 

1.2 

2 

3 

2 

0.  5 

1.0 

0.  5 

0.  75 

1.0 
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Panel  3 


0.  6 

2.1 

0.0  1.5 

2.  3 

2.  3 

0.6  0.0 

2 

1  3 

0.  0 

0.  5 

1.0 

0.  0 

1.0 

GROUP  NO.  3 

Eody  Data 

ITEM. 

_N0Sj _ 7  ZC, 

YC,  AP,  NBE...NT1,  XII, 

RI,  TH1 ,  XIS,  PS 

Eody  1 

o 

o 

0.0 

1.0  1.0 

4 

9  0 

1 

4 

-1.0 

0.0 

1.0  2.0 

3.0 

0.  0 

90.  0 

180.0  270.0 

-2.  0 

-1.3 

-1.4  -1.0 

0.  0 

1.0 

2.  0 

3.  0 

4.0  5.0 

0.  0 

0.  65 

0.95  1.0 

1.0 

1.0 

1.0 

1.0 

0.5  0.0 

Eody  2 

-0.  5 

2.  3 

0.5  1.0 

2 

4  0 

1 

4 

o 

o 

0.  75 

1.5 

45.  C 

1  135.0 

225.0  315.0 

-0.  4 

-0.2 

0.0  0.75 

1.5 

0.  0 

0.  4 

0.5  0.5 

0.  5 

GROUP  NO.  4 

Modal  Spline  Interpolation  Data 

ITEM, 

,N0« _ 1  Control 

Eata 

2 

3  25 

1 

ITEM 

NO, _ 2^4  Control  data,  bodies,  nodes 

Superbody  1 

1 

1 

2  1 

0 

0 

11 

12 

Superbody  2 

1 

2  1 

0 

0 

2 

13 

14 

I  TEM_ 

N2x-2z4  Control  data,  panels,  nodes 

Superpanel  1 

0  4  10 

1 

12  3  4 


Superpanel  2 
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SECTION  VI 


MONITORING  AND  USE  OF  OUTPUT  DATA 

1 .  AERODYNAMIC  DATA  MONITOR  AND  TEST  CASE  OUTPUT 

For  any  two  pre-selected  modes,  e.g.,  symmetric  pitch  and  antisymmetric 
yaw,  a  detailed  set  of  aerodynamic  parameters  are  printed  so  that  the  out¬ 
put  aerodynamics  can  be  properly  monitored.  The  parameters  and  geometric 
data  printed  are  outlined  as  follows: 

(1)  Lifting  pressures  plus  the  locations  of  pressure  points.  (XOC  =  x/c. 
Full  pressures  and  forces  are  printed  for  the  vertical  fin  and 
fuselage  and  not  half  values.) 

(2)  Spanwise  distribution  of  load,  lift  coefficient,  moment  coefficient 
(1/4  chord),  and  center  of  pressure,  versus  spanwise  location. 

(3)  Vertical  and  lateral  running  loads  on  bodies  (force  per  unit  length) 
are  output  versus  longitudinal  location  on  all  bodies.  Output 
given  at  slender  body  element  center  locations. 

(4)  Total  force  and  moment  coefficients  for  the  entire  aircraft  and 
for  lifting  surfaces  only. 

Printed  output  for  the  test  case  described  in  Section  v  is  given  in 
Table  4.  Mode  1  (rigid  body  pitch  of  wing  surface  (not  nacelle)  and  mode 
3  (rigid  body  yaw  of  pylon  and  nacelle)  are  monitored. 

The  following  definitions  are  used  in  the  calculation  of  the  coefficients 
described  above: 

cn  =  c  ^  ACp.AXj  (113) 

j  J 

Cml/4  '  'iSV  <*j  -  >,/«)  <ll4> 

c  J  J 
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Real  C.P. 


Re  Cm]/4/ ReC^  +  X]/4/c 


(115) 


Imag.  C.P.  =  Im  Cm^/ImC^  +  X-^/c 


In  the  equations  to  follow  the  following  constants  are  used: 

A  =  reference  total  area 
c  =  reference  chord 
c  =  local  chord  length 
XM  =  moment  axis  longitudinal  location 
2e  =  spanwise  strip  width 
NB  =  number  of  bodies 
NSTRIP  =  number  of  spanwise  strips 

>r 

£14  =  1/4  chord  point  of  a  box 

J  1  symmetry 

6  =  symmetry  flag  =  » 

(  -1  antisymmetry 

y-  and  z-coordinates  of  the  centerline  of  strip  'j' 


(b) 

^c 

(b) 


=  y-  and  z-coordinates  of  the  centerline  of  slender  body  'I 


M/2  if  y.  =  0  and  cosy.  =  0  and  5/0 
)  J  J 

i  1  otherwise 

(  1/2  if  yc(b)  =  0 
| 1  otherwise 


(116) 
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r>m 


TABLE  4 


TEST  CASE  OUTPUT  LISTING 


11  It  (t 


Kt Tf  r.f r**t am**  $ 


TABLE  4  CONTINUED 


f.iw  P4MIL 


•  ••.FI  V‘,  i  INPUT  V»II|T*  •* 


TABLE  4  CONTINUED 
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2.  TYPICAL  LOAD  DISTRIBUTIONS 


Schematics  of  typical  loadings  are  presented  in  Figure  33  through  36 
as  an  aid  in  the  interpretation  of  the  monitored  output  data.  Specifically 

Figure  33  gives  a  typical  chordwise  loading  for  a  wing  in  steady  pitch 
and  control  surface  deflection.  For  the  pitch  case  a  square  root  singularity 
exists  at  the  leading  edge.  For  the  control  surface  rotation  case  an 
additional  logarithmic  singularity  exists  at  the  control  surface  leading 
edge. 

Figure  34  presents  aspanwise  load  distribution  for;  (a)  a  wing-fuselage 
combination  and  (b)  a  wing-pylon  combination  in  steady  pitch.  The  effect 
of  the  fuselage  is  to  change  the  loading  near  the  wing  fuselage  intersection. 
The  extent  to  which  this  effect  changes  the  wing  load  depends  on  the  body 
diameter  to  wing  span  ratio.  The  pylon  increases  the  inboard  wing  loading. 

The  jump  in  wing  loading,  A,  at  the  pylon  location  is  equal  to  the  pylon  root 
loading. 

Figure  35  presents  a  typical  aerodynamic  center  location  for  a  clean 
wing.  For  swept-back  wings  the  a.c.  moves  aft  near  the  root  and  forward 
near  the  tip.  (The  reverse  is  true  for  wings  with  forward  sweep.)  The  effect 
of  Mach  Number  is  to  move  the  a.c.  further  aft  near  the  root  and  further 
forward  near  the  wing  tip  with  the  cross  over  point  near  the  wing  tip. 

Figure  36  presents  the  vertical  running  load,  load  per  unit  length, 
on  a  wing-fuselage  combination.  This  loading  can  be  thought  of  as  the  sum 
of  two  basic  parts;  a  slender  body  type  of  loading  and  a  lift  carry  through 
type.  The  slender  body  term  changes  most  rapidly  in  regions  where  the  body 
width  or  radius  changes  most  rapidly.  The  wing  lift  carry-over  term  is 
largest  in  the  area  of  the  wing-fuselage  intersection. 


Figure  33  Typical  Chordwise  Loadings 


Figure  34  Typical  Spanwise  Load  Distributions 
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Figyre  35  Typical  Wing  Aerodynamic  Center  Variation 


Figure  36  Typical  Fuselage  Lift  Distribution 


3.  TROUBLE  SHOOTING 


In  this  subsection  a  series  of  problems  are  stated  along  with  their 
possible  cause.  This  list  is  not  all  inclusive  but  covers  many  common 
problem  areas. 

Problem:  Dips  in  span  loading  near  wing-fuselage  intersections  or  in 

the  middle  of  a  wing  or  at  pylon-wing  intersections. 

Possible  Gaps  or  overlap  between  components  or  wing  panels. 

Cause 

Even  a  small  gap  between  wing  and  fuselage  interference  tube 
can  cause  a  significant  dip  in  the  loading  there.  The  same 
holds  true  for  wing-wing  panels  and  wing-pylon  panels. 

Problem:  Irregularities  or  waves  in  the  chordwise  pressure  distribution. 

Possible  A  discontinuous  change  in  box  size  in  the  chordwise  direction. 
Cause 

Or  in  the  case  of  a  fold  line  of  significant  dihedral  change 
(e.g.,  a  wing-pylon  intersection)  the  boxes  do  not  line  up 
in  a  chordwise  direction  on  both  sides  of  the  fold  line.  Or 
there  is  an  overlapping  of  boxes  in  the  chordwise  direction. 

Problem:  Large  or  unreasonable  loads  on  surfaces. 

Possible  Surfaces  may  be  overlapping.  Or  tail  surface  strips  are  not 
Cause 

lined  up  with  wing  strips  in  the  coplanar  or  nearly  coplanar 
case.  Or  intersecting  surfaces  meeting  at  other  than  a  strip 
edge.  In  general  the  problem  may  be  that  control  points  (3/4 
chord  of  box)  lie  too  close  to  singularities;  either  the 
bound  type  (1/4  chord  of  box)  or  the  trailing  type  (lying  on 
strip  edges  starting  at  the  box  and  extending  to  downstream 
infinity) . 
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Problem: 


Possible 

Cause 


Problem: 


Possible 

Cause 


Problem: 


Possible 

Cause 


Problem: 


Possible 

Cause 


The  distribution  has  taken  a  jump  where  there  is  no  pylon 
or  intersecting  surface. 

The  local  chord  length  has  taken  a  similar  jump.  If  intersecting 
surfaces  are  excluded  then  c^c  is  always  continuous.  Thus  if 
c  is  discontinuous  c^  has  the  same  discontinuity. 

Large  or  unreasonable  body  loads. 

The  flow  field  pickup  points  on  the  interference  tube  (given 
by  el.|,  el 2  ...)  lie  at  or  very  near  a  lifting  surface-body 
intersection.  Or  interference  tubes  of  differing  cross-sections 
are  placed  along  the  same  longitudinal  axis. 

Sign  of  some  or  all  of  results  is  incorrect. 

The  panel  may  be  input  upside  down.  That  is,  the  dihedral 
of  the  panel  may  be  different,  by  180°,  from  that  which  was 
intended.  An  example  is  a  vertical  tail  input  from  the  fuselage 
to  tip  has  a  positive  normal  in  negative  y-direction  (dihedral 
=  +90°)  whereas  a  pylon  input  from  wing  down  to  nacelle  has 
a  positive  normal  in  the  positive  y-direction  (dihedral  =  -90°}. 
Or  modes  may  be  input  with  incorrect  signs. 

Body  has  a  loading  even  if  it  is  not  at  angle  of  attack  nor 
intersecting  any  lifting  surface. 

There  is  no  problem.  All  lifting  surface  elements  generate 
some  loading  on  all  bodies  even  though  not  connected  to  that 
body. 
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Problem: 

Possible 

Cause 


Tail  surface  has  unreasonable  loads. 

Tail  surface  lies  within  an  interference  tube.  Only  one 
interference  tube  is  allowed  on  any  one  longitudinal  axis, 
e.g.,  two  interference  tubes,  one  for  the  wing  and  one  for 
the  tail,  can  not  be  used  for  a  single  fuselage.  Or  tail  strips 
do  not  line  up  with  wing  strips  in  the  coplanar  or  nearly 
coplanar  case. 
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SECTION  VII 


SUBROUTINE  DESCRIPTION 


1.  INTRODUCTION 

A  large  portion  of  the  computer  program  being  introduced  into  the 
VI BRA-6  system  is  documented  in  Reference  4  (Vol.  II  of  Part  II),  and 
need  not  be  repeated  here.  However  new  subroutines  have  been  developed, 
old  ones  modified  and  some  have  had  their  argument  lists  changed  to  eliminate 
common  blocks  of  data.  Twenty-nine  of  the  forty-six  subroutines  of  the 
program  have  had  some  sort  of  change,  and  fourteen  subroutines  have  been 
added. 

No  limits  on  dimension  sizes  have  been  assigned  since  a  dynamic  core 
feature  is  to  be  used  for  the  entire  VIBRA-6  system.  The  working  array 
size  will  depend  on  the  problem  solved.  Formulas  dealing  with  the  working 
array  size  are  furnished  in  Volume  I. 

The  listing  of  the  program  and  subroutines  is  given  in  Volume  III. 
Descriptions  of  each  of  the  thirty-four  modified  and  changed  subroutines 
follow  in  alphabetic  order. 


102 


2.  DETAILED  DESCRIPTION  OF  SUBROUTINES 

2.1  SUBROUTINE  AERO  (MFIX1,  MFIX2,  NMSYM,  NTOT,  NOUT,  NW,  NEWBFZ 


NEWBFY ,  NSTOT ,  CBSPAN,  DCP ,  FZ,  FY,  CN,  CM,  SPLD, 
ISSTR ,  NSBEA,  NBARAY ,  NCARAY ,  YB,  ZB,  XIS1  .  XIS2, 
CG,  CS,  EE,  SG.  YS,  ZS,  XIC,  XIJ,  DELX ,  COORD,  CZB. 
CYB,  CNB ,  CMB ,  CPR,  CPI) 


Functional  Descirption 

Subroutine  AERO  computes  and  prints  the  pressures  and  the  aerodynamic 
derivatives  for  two  fixed  modes,  MFIX1  and  MFIX2.  It  is  called  once 
for  each  reduced  frequency  from  subroutine  SDLM. 


Input  Output  Variables 


MNEMONIC 

MFIX1 

MFIX2 

NMSYM 

NT0T 

N0UT 

NW 

NEWBFZ 

NEWBFY 

NSTOT 

CBSPAN 


DESCRIPTION 

MFIX1  and  MFIX2  are  two  mode  numbers  for  which 

aerodynamic  parameters  are  calculated 

The  total  number  of  symmetric  modes 

The  total  number  of  unknowns 

Data  set  number  of  the  system  output  data  set 

Tape  number  containing  the  pressures  for  all 

points  in  all  modes 

Tape  number  containing  the  body  element 
z-forces 

Tape  number  containing  the  body  element  y-forces 
Maximum  dimension  of  body  force  array  or  sectional 
lift  coefficient  array 

Reference  length  used  in  the  calculation  of  span 
loads 
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MNEMONIC 

DCP 

FZ 

FY 

CN 

CM 

SPLD 

CZB 

CYB 

CNB 

CMB 

CPR 

CPI 


DESCRIPTION 


Complex  array  of  pressures 
in  one  mode 

Body  element  z-forces  1 
Body  element  y-forces  ; 
Lift  coefficients  ) 


for  all  points  (boxes) 


in  one  mode 


for  all  strips  in  one  mode 


Moment  coefficients  ' 

Span  loads 

Total  lift  coefficients  for  bodies,  z-  and  y- 

directions 

Total  yawing  Y 

>  moment  coefficients  for  bodies 
Total  pitching  ) 

Centers  of  pressure  (real  -  and  imaginary  parts) 
for  lifting  surface  strips 


The  following  variables  are  described  under  Basic  Variables  in 
Section  viI-3: 

ISSTR,  NSBEA,  NBARAY,  NCARAY ,  YB,  ZB,  XIS1,  XIS2,  CG,  CS,  EE,  SG, 
YS,  ZS,  XIC,  XI J ,  DELX,  COORD. 


Calling  Subroutine  SDLM 
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2.2  SUBROUTINE  BFSMAT ( LENGTH .  NE,  NB.  NP,  NTP ,  NTOTAL,  IQ,  IOA,  FMACH, 

CBAR,  KR,  YB,  ZB,  YS,  ZS,  X,  DELX,  EE,  XIC,  SG,  CG,  AR,  RIA, 

AO,  XIS1,  X I S2 ,  AVR,  NSARAY,  NCARAY,  NSBEA,  NBEA1 ,  NAS,  NASB, 
BPS,  BFSA,  SCALER) 

Functional  Description 

This  subroutine  is  the  basic  calling  routine  that  computes  the  [FZ]^  and 
[FY]^)  matrices  for  both  symmetric  and  antisymmetric  conditions.  These 
matrices  are  written  on  two  tapes,  in  row  order,  alternating  first  a 
row  of  FZ,  then  a  row  of  FY  for  symmetry  is  written  on  tape  (scratch  unit) 
10,  and  similarly,  a  row  of  FZ,  then  a  row  of  FY  for  antisymmetry  is 
written  on  another  tape,  IOA. 


Input  Output  Variables 


MNEMONIC 

LENGTH 

NE 

NB 

NP 

NTP 

NT0TAL 

10 

I0A 

FMACH 

CBAR 

KR 

YB 


DESCRIPTION 

The  total  number  of  unknowns  +  the  total  number 

of  z-  and  y-oriented  slender  body  elements 

Ground  effect  flag 

Number  of  bodies 

Number  of  panels 

Number  of  lifting  surface  boxes 

NTP  +  Total  number  of  z-  and  y-oriented  interference 
body  elements 

Tape  number  on  which  the  BFS-matrix  is  written 

Tape  number  on  which  the  BFSA-matrix  is  written 

Mach  Number 

Reference  chord 

Reduced  frequency 

Array  of  y-coordinates  of  bodies 


Array  of  z-coordinates  of  bodies 
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ZB 


MNEMONIC 


DESCRIPTION 


YS 

ZS 

X 

DELX 

EE 

XIC 

SG 

CG 

AR 

RIA 

AO 

XIS1 

XIS2 

AVR 

NSARAY 

NCARAY 

NSBEA 

NBEA1 

NBEA2 

BFS 

BFSA 


Array  of  y-coordinates  of  lifting  surface  strips 
Array  of  z-coordinates  of  lifting  surface  strips 
Array  of  3/4-chord  locations  of  boxes  and  1/2 -chord 
for  body  elements 

Array  of  lengths  of  panel  boxes  and  body  elements 
Array  of  the  semiwidths  of  strips 
Array  of  1/4-chord  locations  of  all  boxes 
Array  of  the  sines  of  the  dihedral  angles  of  strips 
Array  of  the  cosines  of  the  dihedral  angles  of  strips 
Array  of  the  cross-sectional  aspect  ratics  of  the  bodies 
Array  of  the  radii  of  interference  body  elements 
Array  of  the  slender  body  element  radii 

Leading  Y  e(jge  ]ocafj0n  0f  the  slender 

Trailing  I  body  elements 

Array  of  the  average  radii  of  all  bodies 

Array  of  the  number  of  strips  per  panel 

Array  of  the  number  of  chordwise  boxes  per  panel 

Array  of  the  number  of  slender  body  elements  per  body 

Array  of  the  number  of  interference  body  elements 

per  body 

z-y  orientation  flag  array  per  bod,'  fixed  as  2) 

One  row  of  the  z-  and  y-forces  in  the  BFS  matrix 
(symmetry) 

One  row  of  the  z-  and  y-forces  in  the  BFSA  matrix 
(antisymmetry) 

Nondimensional  image  radius  (input) 


SCALER 


SDLM 


Calling  Subroutine 
Called  Subroutine 


FWMW 
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2.3  SUBROUTINE  CTLS  (NEAASS,  NDAS,  ND,  ELXI,  ELY I ,  ELZI .  XA.  YA,  ZA 


XB,  YB,  ZB,  XAR,  YAR,  ZAR,  XBR,  YBR,  ZBR,  XI,  ETA,  NTP6) 
Functional  Description 

This  subroutine  computes  the  ?,  and  n-coordi nates  of  a  generated  node 
on  the  elastic  axis  of  a  control  surface,  for  which  no  input-node  exists, 
from  an  associated  node  (NDAS)  of  the  associated  elastic  axis  (NEAASS). 
(See  Figure  14). 


Input  Output  Variables 


MNEMONIC 

NEAASS 

NDAS 

ND 


ELXI 

ELYI 

ELZI 

XA, 

YA,  ZA  | 

GO 

X 

YB,  ZB  ' 

XAR, 

YAR,  ZAR 

XBR, 

YBR,  ZBR 

DESCRIPTION 

Sequence  number  of  the  associated  elastic  axis 
(i.e.,  superpanel) 

Associated  node  number 

Generated  node  number  for  the  current  control 

surface  elastic  axis;  numbering  begins  with  ND  =  1 

for  each  control  surface,  (1<ND s  NODASS,  where 

NODASS  is  card  input) 

Arrays  of  the  x-,  y-,  z-coordinates 

of  all  input  nodes  on  all  elastic 

axes  /  \ 

‘inboard  ( 

Coordinates  of  the  |outboardj  edge  of  elastic 

axis  for  control  surface  ,  j 

‘inboard  (, 

Coordinate  arrays  of  the  joutboardj  edges  of 
all  input  elastic  axes 


I  , 


MNEMONIC 
XI,  ETA 
NTP6 

Calling  Subroutine 


DESCRIPTION 

£  ,  n  coordinate  arrays  of  the  generated  nodes 
Reference  number  of  the  system  output  unit 


NEWH 
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2.4  SUBROUTINE  DOTP  (NTPH,  NTPDH ,  NTPH4,  MASTH,  MASTDH,  MASTH4 ,  IFNEWH, 

NSYM,  NASYM,  NODES,  MODES,  NB,  I  ROW,  PHIN,  PHIZ,  PHIY,  COL, 
WORK) 

Functional  Description 

This  subroutine  organizes  the  input/output  arguments  to  subroutine  ORGN 
for  the  computation  steps  of  the  [H]  matrices  for  both  symmetric  and 
antisymmetric  conditions  for  lifting  surface  panels  as  well  as  bodies. 

Input  Output  Variables 

MNEMONIC  DESCRIPTION 


NTPH 

Tape  number  containing  the  h  arrays;  output  of 

subroutine  SPLINE 

NTPDH 

Tape  number  containing  the  dh/dx  arrays;  output 

of  subroutine  SPLINE 

NTPH4 

Tape  number  containing  the  h^4c  arrays;  output 

of  subroutine  SPLINE 

MASTH 

1  Output  tapes  containing  the  results  of  the  dot 

MASTDH 

(  _ 

/  products,  (hi,  {dh/dx}  and  {hy4c> 

MASTH4 

'  respectively 

IFNEWH 

Flag  to  select  degrees  of  freedom;  see 

subroutine  SDLM 

NSYM 

The  number  of  symmetric  modes 

NASYM 

The  number  of  antisymmetric  modes 

NODES 

The  number  of  nodes 

MODES 

The  number  of  modes 

MB 

The  number  of  bodies 
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MNEMONIC 

I  ROW 

PHIN 

PHIZ 

PHI  Y 

COL 

WORK 


DESCRIPTION 

Maximum  number  of  row  elements  in  the  work  array 
Modal  deflection  matrix  for  panels 
Modal  deflection  matrix  for  bodies  in  z-direction 
Modal  deflection  matrix  for  bodies  in  y-direction 
One  column  of  the  h,  dh/dx  or  hy^c  matrices 
Output  of  subroutine  ORGN  containing  either  of  the 
matrices  [h],  [dh/dx]  or  [h-j^]  which  are  saved 
on  the  three  tapes  MASTH,  MASTDH  and  MASTH4 
respectively 


Calling  Subroutine  SDLM 


Called  Subroutine 


ORGN 


2.5  SUBROUTINE  DPPS  (K,  KS,  I,  Jl,  J2,  SGR.  CGR,  REFC,  FMACH,  YS.  ZS, 
NBARAY ,  NCARAY,  DT,  DTA,  YB,  ZB,  ARB,  AVR,  XLE,  XTE , 

X,  CG,  EE,  SG,  XIC,  DELX ,  XLAM,  SCALFR) 

Functional  Description 

Subroutine  DPPS  prepares  the  variables  necessary  for  the  computation 
of  one  row  of  the  DPP  submatrix  (for  symmetry)  and  one  row  of  the 
DPPA  submatrix  (for  antisymmetry) .  The  resulting  two  matrix  rows  are 
returned  to  subroutine  GEND  via  the  argument  list  in  arrays  DT  and 
DTA  respectively. 


Input  Output  Variables 


MNEMONIC 

K 

KS 

I 

Jl,  J2 

SGR,  CGR 

DT 

OTA 

SCALER 


DESCRIPTION 

Panel  number  in  which  the  receiving  point 
1 i *  lies 

Strip  number  containing  1  i ' 

Receiving  point  index 

DO  loop  delimiters  for  the  number  of  elements 

in  one  row  of  the  submatrices 

Sine  and  cosine  of  the  receiving  strip  dihedral 

angle 

A  row  of  the  DPP  submatrix 
A  row  of  the  DPPA  submatrix 
Nondimemsional  image  radius 
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The  following  variables  are  described  under  Basic  Variables  in 
Section  VII ,  Part  3: 

REFC,  FMACH,  YS,  ZS,  NBARAY ,  NCARAY ,  YB,  ZB,  ARB,  AVR,  XLE,  XTE,  X, 
CG,  EE,  SG,  XIC,  DELX,  XLAM 

Calling  Subroutine  GEND 

Called  Subroutine  SUBP 
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5  SUBROUTINE  DPZY(KB,  IZ,  I,  J1,  J2,  IFIRST,  HAST,  REFC,  FMACH, 

YB,  ZB,  AVR,  ARB,  TH1A ,  TH2A,  NT12,  NBARAY ,  NCARAY , 

NZYKB,  DPZ ,  DPY ,  DPZA,  DPYA,  XLE,  XTE,  X,  CG,  EE, 

SG,  ys,  zs,  xic,  delx,  xlam,  scaler 

Functional  Description 
Subroutine  DPZY  prepares  the  variables  necessary  for  the  computation 
of  one  row  of  either  the  DPZ  or  DPY  submatrix  (for  symmetry),  and 
one  row  of  either  the  DPZA  or  DPYA  submatrix  (for  antisymmetry)  and 
calls  subroutine  SUBP  in  a  double  DO  loop  for  each  element  of  the 
row  to  perform  the  necessary  summation  around  an  interference  body 
segment.  The  resulting  matrix  rows  are  returned  to  subroutine  GEND 
via  the  argument  list. 


Input  Output  Variables 

MNEMONIC  DESCRIPTION 

KB  Body  number  in  which  the  receiving  point 

' i '  lies 

IZ  Body  element  number  of  body  KB  in  which 

’ i  ’  lies 


I 

Jl,  J2 

IFIRST 

HAST 

NZYKB 

SCALER 


Receiving  point  index 

DO  loop  delimiters  for  the  number  of  elements 
in  a  row  of  the  submatrices 

Sequence  number  of  the  first  element  in  body  KB 
Sequence  number  of  the  last  element  in  body  KB 
z-y  flag  of  body  KB 
Nondimensional  image  radius 
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The  following  variables  are  described  under  Basic  Variables  in 
Section  VI I -3 : 


REFC,  FMACH,  YB,  ZB, 
DPY ,  DPZA,  DPYA,  XLE 

Calling  Subroutine 
Called  Subroutine 


AVR,  ARB,  TH1A,  NT12,  NBARAY ,  NCARAY ,  DPZ, 
XTE,  X,  CG,  EE,  SG,  YS,  ZD,  XIC,  DELX,  XLAM. 

GEND 

SUBP 
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2.7  SUBROUTINE  DYPZ(KB,  KS,  LS,  IZ,  I.  Jl,  J2,  NYFLA6,  FLND,  FLNE, 

SGR,  CGR,  REFC,  FMACH,  KR,  ARB,  NBEA,  LBO,  LSO,  JBO, 

DT,  DTA,  YB,  ZB,  RIA,  X.  YS.  ZS,  DELX) 

Functional  Description 

Subroutine  DYPZ  prepares  the  variables  necessary  for  the  computation 
of  one  row  of  either  the  DYP,  or  the  DYZ,  or  the  DYY  submatrix  (for 
symmetry)  and  one  row  of  either  the  DYPA,  or  the  DYZA,  or  the  DYYA 
submatrix  (for  antisymmetry).  In  either  case  it  calls  subroutine 
SUBB  in  a  DO  loop  for  each  element  of  a  row;  latter  is  returned  to 
subroutine  GEND  via  the  argument  list. 


Input  Output  Variables 


MNEMONIC 

KB 

KS 

LS 

IZ 

I 

Jl,  J2 


DESCRIPTION 

Body  number  in  which  the  receiving  point 
' i ’  lies 

Index  of  the  receiving  point  y,  z  coordinates 
Index  of  the  sending  point  y,  z  coordinates 
Body  element  number  in  which  1 i '  lies 
Receiving  point  index 

DO  loop  delimiters  for  the  number  of  elements 
in,  a  row  of  the  submatrices 


NYFLAG 


NYFLAG=0  for  DYP,  DYZ,  DYPA  and  DYZA 
NYFLAG=1  for  DYY  and  DYYA 


FLND 

FLNE 

SGR,  CGR 
LBO 


Not  used 

e,  ground  effect  flag 

Sine  and  cosine  of  receiving  point  dihedral  angle 
Sequence  number  of  first  body  with  y  orientation 
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MNEMONIC 


DESCRIPTION 


LSO  Index  for  the  y,  z  coordinates  of  the  first 

y-oriented  body 

JBO  Sending  point  index  for  first  y-oriented  body 

element 

DT  A  row  of  either  of  the  submatrices  DYP,  DYZ 

or  DYY 

DTA  A  row  of  either  of  the  submatrices  DYPA,  DYZA 

or  DYYA 

The  following  v  riables  are  described  under  Basic  Variables  in 

Section  VII-3: 

REFC,  FMACH,  KR,  ARB,  NBEA,  YB,  ZB,  RIA,  X,  YS,  ZS,  DELX 

Calling  Subroutine  GEND 

Called  Subroutine  SUBB 
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2.8  SUBROUTINE  DZPY(KB,  KS,  LS,  IZ,  I,  Jl,  J2,  NYFLAG,  FLND,  FLNE, 

SGR,  CGR,  REFC ,  FMACH ,  KR.  ARB,  NBEA,  DT,  DTA,  YB, 

ZB,  RIA,  X,  YS,  ZS,  DELX) 

Functional  Description 

Subroutine  DZPY  prepares  the  variables  necessary  for  the  computation 
of  one  row  of  either  the  DZP,  or  the  DZZ  or  the  DZY  submatrix  (for 
symmetry),  and  one  row  of  either  the  DZPA,  or  the  DZZA  or  the  DZYA 
submatrix  (for  antisymmetry).  In  either  case  it  calls  subroutine 
SUBB  in  a  double  DO  loop  for  each  element  of  a  row;  latter  is 
returned  to  subroutine  GEND  via  the  argument  list. 

Input  Output  Variables 

MNEMONIC  DESCRIPTION 

KB,  KS 
LS 
17. 

I  See  Description  given 

Jl ,  J2  in  Section  2.6 

NYFLAG 

FLND 

FLNE 

SGR,  CGR 

DT  A  row  of  either  of  the  submatrices 

DZP,  or  DZZ,  or  DZY 

DTA  A  row  of  either  of  the  submatrices 

DZPA,  or  DZZA,  or  DZYA 
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The  following  variables  are  described  under  Basic  Variabl 
Section  VII-3: 

REFC,  FMACH ,  KR,  ARB,  NBEA,  YB,  ZB,  RIA,  X,  YS,  ZS,  DELX. 


Calling  Subroutine 
Called  Subroutine 


GEND 

SUBB 


2.9  SUBROUTINE  DZYMAT(D,  DA,  NFB,  NLB,  NTZYS,  IDZDY ,  NTAPE,  NTAP2, 

XP,  BETA,  IPRNT,  NBEA,  NS BE,  NC.  NS,  AO,  YB,  ZB,  AR. 

XIS1,  XIS2 ,  CG,  SG,  YP,  ZP) 

Functional  Description 

This  subroutine  sets  up  the  argument  list  for  the  calculation  of  each 
row  of  the  Dz  and  matrices  both  in  symmetry  and  antisymmetry  and 
then  calls  subroutine  ROWDYZ. 


Input  Output  Variables 


MNEMONIC 

D 

DA 

NFB 

NLB 

NTZYS 

IDZDY 

NT  APE 

NTAP2 

XP 

BETA 

IPRNT 


DESCRIPTION 

Complex  array  containing  a  row  of  the  Dz  or  D^ 
matrices  in  symmetry 

Complex  array  containing  a  row  of  the  Dz  or  Dy 

matrices  in  antisymmetry 

Number  of  the  first  body 

Number  of  the  last  body 

Number  of  slender  body  elements 

Flag  indicating  whether  the  Dz,  or  the  D^  matrix 

is  calculated 

Tape  number  (scratch  unit)  on  which  the  output 
matrix  rows,  D,  are  saved 

Tape  number  (scratch  unit)  on  which  the  output 
matrix  rows,  DA,  are  saved 

x-coordinate  of  lifting  surface  box  control  points 
6  =  y 1  -  M2 ,  M  =  Mach  Number 
Not  used 
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The  following  variables  are  described  under  Basic  Variables  in 
Section  VII-3: 

NBEA,  NSBE,  NC,  NS,  AO,  YB,  ZB,  AR,  XIS1  ,  XIS2,  CG,  SG,  YP,  ZP. 


Calling  Subroutine 
Called  Subroutine 


SB 

ROWDYZ 
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2.10  SUBROUTINE  FVIMW(NDX,  NE,  SGS,  CGS,  IRB,  AO,  ARB,  XBLE,  XBTE,  YB,  ZB, 
XS.  YS,  ZS,  NAS,  NASB,  KR.  BETA2 ,  CBAR,  AVR,  FWZ,  FUZA,  FWY , 
FWYA,  IF1 ,  IPRNT ,  SCALER) 

Functional  Description 
Given  a  unit  pressure  doublet, this  subroutine  calculates  the  effect  of 
this  doublet  plus  any  contributions  due  to  images,  symmetry  plane  and 
ground  effect  on  a  given  body  both  in  symmetric  and  antisymmetric 
conditions. 


Input  Output  Variables 
MNEMONIC 


NDX 

NE 

SGS 

CGS 

IRB 

AO 

ARB 

XBLE 

XBTE 

YB 

ZB 

XS 

YS 

ZS 


DESCRIPTION 


Not  used 

Ground  effect  flag 
Sine 

•  of  sending  point  dihedral  angle 

Cosine 

Number  of  the  receiving  body 
Radius  of  body  element 
Cross-sectional  aspect  ratio  of  body 
Leading 
Trailing 

Array  of  y-coordinates  of  bodies 
Array  of  z-coordinates  of  bodies 
x-coordinate 


f  edge  location  of  slender  body  element 


y-coordinate 

z-coordinate 


L  of  sending  point 
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MNEMONIC 


DESCRIPTION 


NAS 

Number  of  bodies 

NASB 

Array  of  body  numbers 

KR 

Reduced  frequency 

BETA2 

B2  =  1-  M2,  M  =  Mach  Number 

CBAR 

Reference  chord 

AVR 

Array  of  the  average  radii  of  bodii 

FWZ 

Output  z-forces  in  symmetry 

FWZA 

Output  z-forces  in  antisymmetry 

FWY 

Output  y-forces  in  symmetry 

FWYA 

Output  y-forces  in  antisymmetry 

IF1 

Orientation  flag  of  receiving  body 

IPRNT 

Not  used 

SCALER 

Nondimensional  image  radius 

Calling  Subroutine  BFSMAT 

Called  Subroutine  FZY2 ,  SUBI 


2.11  SUBROUTINE  GEND(NPRINT,  NTAPE .  NTOT,  IQPT,  MASTDT,  DT,  DPZ,  DPY , 
DTA,  DPZA.  DPYA.  I  FLA,  NBEA,  NT!  2,  NBARAY ,  NCARAY ,  YB 


ZB,  ARB,  AVR,  RIA,  XLE,  XTE,  TH1A,  TH2A ,  X,  CG,  EE, 

SG,  YS,  ZS,  XIC,  YIN,  ZIN,  DELX,  XLAM.  SCALER 1 

Functional  Description 

Subroutine  GEND  generates  the  nine  submatrices  of  the  symmetric  downwash 
factor  matrix,  [DT],  and  the  nine  submatrices  of  the  antisymmetric 
downwash  factor  matrix,  [DTA].  These  two  matrices  are  written  on  a 
utility  tape  in  row  order,  one  row  of  the  DT  matrix  and  one  row  of 
the  DTA  matrix  as  one  record. 


Input  Output  Variables 


MNEMONIC 

NPRINT 

NTAPE 

NTOT 

IOPT 


MASTDT 

DT 

DPZ 

DPY 

DTA 


DESCRIPTION 

Print  flag  for  DT  and  DTA  matrix  rows; 

NPRINT  =  0  for  no  print 

Array  of  utility  tape  numbers 

Total  number  of  unknowns 

Option  flag  for  saving  DT  and  DTA  on  output 

tape  for  future  use,  when  IOPT  =  1; 

IOPT  =  0  otherwise 

Output  tape  number  on  which  the  DT  and  DTA 
rows  are  saved  when  IOPT  =  1 
Array  containing  one  row  of  the  DT  matrix 
(symmetric  downwash  factors) 

Array  containing  one  row  of  the  DPZ  submatrix 

Array  containing  one  row  of  the  DPY  submatrix 

Array  containing  one  row  of  the  DTA  matrix 

(antisymmetric  downwash  factors) 
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MNEMONIC 


DESCRIPTION 


Array  containing  one  row  of  the  DPZA  submatrix 
Array  containing  one  row  of  the  DPYA  submatrix 
Nondimensional  image  radius 

The  following  variables  are  described  under  Basic  Variables  in 
Section  VI 1-3: 

IFLA,  NBEA,  NT1 2,  NBARAY,  NCARAY ,  YB,  ZB,  ARB,  AVR,  RIA,  XLE,  XTE, 
TH1A,  X,  CO,  EE,  SG,  YS,  ZS,  XIC,  YIN,  ZIM,  DELX ,  XL AM. 

Calling  Subroutine  SDLM 

Called  Subroutines  DPPS,  OPZY,  DZYP,  DYPZ 


DPZA 

DPYA 

SCALER 
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2.12  SUBROUTINE  GENF  (NMODES,  NTQT,  NTP,  NB,  YB,  NSBETO.  IPRINT,  LINES, 

NTP6,  NUTL1  ,  NUTL2,  NTPH,  NTPH4,  NTFORC,  NTGF,  MASTCP, 

MASTFZ ,  MASTFY ,  JKR,  H,  FMULT,  NASYM,  DELA,  XIS1  ,  XIS2, 

NSBEA  FORCE,  GF,  WORK,  NG) 

FUNCTIONAL  DESCRIPTION 

This  subroutine  calculates  element  forces  and  generalized  forces  for  all 
modes  and  gusts.  The  pressures  ACp  are  multiplied  by  the  box  area  aA, 
to  obtain  the  forces 

[D0  :  Fg]  =  aA  [ACp] 

D0  represents  the  forces  due  to  modes  and  Fg  represents  the  forces  due 
to  gusts.  This  operation  is  done  for  lifting  surfaces  (panels),  and 
bodies  and  for  both  symmetric  and  antisymmetric  modes  and  gusts. 

The  generalized  forces  are  calculated  using  these  forces: 

^  :<?6]  -  [h1/4J  t  cm  :  fg] 

Again  this  operation  is  repeated  for  lifting  surfaces  and  bodies  and 
symmetric  and  antisymmetric  modes.  The  forces  and  generalized  forces 
are  written  on  tape. 


Input  Output  Variables 


MNEMONIC 

NM0DES 

NT0T 

NTP 

NB 

YB 


DESCRIPTION 
Array  of  mode  numbers 
Maximum  dimension  of  the  force  array 
Number  of  lifting  surface  boxes 
Number  of  bodies 

Array  of  the  y-coordinates  of  body  centerlines 
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MNEMONIC 


DESCRIPTION 


NSBETO 

IPRINT 

LINES 

NTP6 

NUTL1 

NUTL2 

NTPH 

NTPH4 

NTFORC 


NT6F 

MASTCP 

MASTFZ 

MASTFY 

JKR 

H 

FMULT 

NASYM 

DELA 

XIS1 

XIS2 

NSBEA 

FORCE 

GF 


Total  number  of  slender  body  elements 
Print  flag 

Maximum  lines  of  print  per  page 

Data  set  number  of  the  system  output  data  set 

Utility  (scratch)  tape  number 

Utility  (scratch)  tape  number 

Tape  number  containing  the  [h]  matrix 

Tape  number  containing  the  [hj^]  matrix 

Tape  number  containing  the  output  items  of  the 

aero  module  in  the  VIBRA-6  program  (includes 

point  forces  and  generalized  forces) 

Tape  number  containing  the  generalized  forces 
Tape  of  pressures  (on  lifting  surface  boxes) 

Tape  of  z-body  element  forces 
Tape  of  y-body  element  forces 
Index  of  current  reduced  frequency 
Array  of  deflections 
Array  for  temporary  use 
Number  of  antisymmetric  modes 
Array  of  lifting  surface  areas 
Leading  i 

,  edge  location  of  slender  body  elements 
Trailing  ) 

Array  of  the  number  of  slender  body  elements  per 
body 

One  column  of  the  point  force  matrix 
One  column  of  the  generalized  force  matrix 
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MNEMONIC 


DESCRIPTION 


W0RK  Complex  work  array 


Calling  Subroutine  SDLM 
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SUBROUTINE  GEOM  (DELA,  I FLA.  NBEA.  NT12,  NAS.  NASB.  ISSTR.  NSSTR. 

NSBEA,  NBARAY ,  NCARAY ,  NSARAY,  AO,  YB,  ZB,  ARB.  AVR. 

AOP ,  RIA,  XLE,  XTE ,  TH1A,  XIS1 ,  XIS2,  X.  CS.  CS,  EE, 

SG,  YS,  ZS,  XIC ,  XU.  YIN.  ZIN,  DELX,  XL  AM,  COORD, 
NTAERO,  TH,  TAU,  XII.  RI ,  XIS.  RS,  GMA,  DYS,  DZS,  GMAR, 
XII,  XI2,  ETA1 ,  ETA2,  ZETA1 ,  ZETA2 ,  XC,  ETA,  ZETA, 

ETAS,  ZETAS,  I HD,  TEMP,  I TEMP) 


Functional  Description 

Subroutine  GEOM  calculates  the  geometry  array  elements  for  all  lifting 
surface  panels  and  all  interference  body  and  slender  body  elements. 

The  basic  geometry  arrays  are  saved  in  core  and  are  used  throughout 
the  program.  Those  coordinate  arrays  that  are  required  subsequent  to 
the  aero  module  of  VIBRA-6,  are  saved  on  the  output  tape  NTAERO. 


Input  Output  Variables 


MNEMONIC 

NTAERO 

TH 

TAU 

XII 

RI 

XIS 

RS 


DESCRIPTION 

Output  tape  number  containing  the  coordinate 
arrays  for  the  case 

Fractional  chordwise  divisions  for  panel 
Fractional  spanwise  divisions  for  panel 
X-coordinates  of  the  interference  body 
element  endpoints 

Radii  of  the  interference  body  element  endpoints 
X-coordinates  of  the  slender  body  element 
endpoints 

Radii  of  the  slender  body  element  endpoints 
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MNEMONIC 


DESCRIPTION 


GMA 

DYS,  DZS 
GMAR 
XII,  XI2 

ETA1 ,  ETA2 


ZETA1 ,  ZETA2 


XC 

ETA 

ZETA 

ETAS 

ZETAS 

IHD 

TEMP 

I  TEMP 


Dihedral  angles  for  all  panels 
Ay,  az  of  panel  strips 
Array  of  dihedral  angles  for  all  strips 
Inboard,  outboard  x-coordinates  of  the  1/4-chord 
location  of  all  lifting  surface  boxes 
Inboard,  outboard  y-coordinates  of  all  lifting 
Surface  boxes;  also,  y-coordinates  of  all 
slender  body  elements  (centerline) 

Inboard,  outboard  z-coordinates  of  all  lifting 
surface  boxes;  also,  z-coordinates  of  all  slender 
body  elements  (centerline) 

X-coordinates  of  lifting  surface  box  corners 
Y-coordinates  of  lifting  surface  box  centerlines 
Z-coordinates  of  lifting  surface  box  centerlines 
Y-coordinates  of  slender  body  element  centerlines 
Z-coordinates  of  slender  body  element  centerlines 

Temporary  storage  for  selected  constants  and 
geometry  arrays  of  a  case 


The  following  variables  are  described  under  Basic  Variables  in 
Section  VII-3: 

DELA,  I  FLA,  NBEA,  NT12,  NAS,  NASB,  ISSTR,  NSSTR,  NSBEA,  NBARAY ,  NCARAY 
NSARAY,  AO,  YB,  ZB,  ARB,  AVR,  AOP,  RIA,  XLE,  XTE,  TH1A,  XIS1  ,  XIS2,  X, 
CG,  CS,  EE,  SG,  YS,  ZS,  XIC,  XIJ,  YIN,  ZIN,  DELX,  XLAM,  COORD. 

Calling  Subroutine  SDLM 

Called Subroutine  GRUP 
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2.14  SUBROUTINE  GRUP(NS,  J2,  IOLD,  IGRUP,  MAXSSTR,  ISSTR,  NSSTR,  COORD) 
Functional  Description 

This  subroutine  calculates  the  superstrip  number  of  each  strip  on  all 
lifting  surface  panels  and  saves  these  in  an  array  to  be  used  in  sub¬ 
routine  AER0  for  the  span  load  calculations.  Input/output  to  this 
subroutine  is  via  the  argument  list. 

Input  Output  Variables 

MNEMONIC  DESCRIPTION 


NS 

Number  of  strips  in  panel 

J2 

Strip  number  of  last  strip  in  panel 

IOLD 

Group  number  of  preceding  panel; 

I0LD=O  for  the  first  panel 

IGRUP 

Group  number  of  current  panel 

MAXSTR 

Number  of  superstrips 

ISSTR 

Superstrip  number  of  each  strip 

NSSTR 

Number  of  strips  per  superstrip 

COORD 

Spanwise  coordinate  of  strips 

Calling  Subroutine 

GEOM 
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2.15  SUBROUTINE  GUST  (NWT,  NWTA,  NTOT,  NBQX,  NB,  NBARAY ,  NCARAY ,  NBEA, 

FMACH,  KR,  CBAR,  X,  YS,  ZS,  SG,  CG,  CR,  NG,  IPRINT , 

LINES,  NOUT,  MGS ,  WGA,  COL,  NDU,  NMT ,  NASYH) 

Functional  Description 

This  subroutine  computes  the  gust  boundary  conditions  for  thirteen  symmetric, 
and  thirteen  antisymmetric  modes  according  to  the  equations  described  in 
Section  II.  (See  Figure  5). 


Input  Output  Variables 
MNEMONIC 


NCARAY 

NBEA 

FMACH 


DESCRIPTION 


NWT 

Tape  number  on  which  the  symmetric  gust  boundary 

conditions  are  saved 

NWTA 

Tape  number  on  which  the  antisymmetric  gust 

boundary  conditions  are  saved 

NTOT 

Number  of  unknowns 

NBOX 

Number  of  boxes  (elements)  on  lifting  surfaces 

NB 

Number  of  bodies 

NBARAY 

nbap  =  ^  (nc^  (ns..),  where 

in 

panel 

i 


i  =  l 

p  =  panel  number 

nc.j  =  number  of  chordwise  boxes  1 
ns.j  =  number  of  spanwise  strips  J 

nc  =  number  of  chordwise  boxes  per  panel 
Number  of  interfence  body  elements  per  body 
Mach  Number 
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MNEMONIC 


DESCRIPTION 


KR 

CBAR 

X 

YS 

ZS 

SG 

CG 

CR 

NG 

I PR I NT 

LINES 

NOUT 

WGS 

WGA 

COL 

NDW 

NMT 

NASYM 


Reduced  frequency 
Reference  chord 

3/4-chord  x-coordi nates  of  all  lifting  surface  boxes 
and  1 /2-chord  x-coordl nates  of  body  elements 
y- coordi nates  1  of  strip  centerlines  and 
z-coordi nates  )  of  body  axes 
sine  1 

\  of  dihedral  angles  of  strips 
cos i ne ) 

Array  of  the  direction  cosines  for  the  gust 
boundary  conditions 
Number  of  gusts 
Print  flag 

Maximum  lines  of  print  per  page 

Data  set  number  of  the  system  output  data  set 

One  column  of  the  symmetric  gust  boundary  conditions 

One  column  of  the  antisymmetric  gust  boundary 

conditions 

Array  for  temporary  use 

Tape  number  containing  the  Incremental  normalwash 
aW  due  to  slender  bodies 

Total  number  of  symmetric  and  antisymmetric  modes 
Number  of  anti syrmetri c  modes 


Calling  Subroutine  SDLM 
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2.17 


SUBROUTINE  MATMUL(NW,  NPSTAP,  NBFM,  NEWBFZ,  NEWBFY,  NWORK ,  NTSBE, 

NSBETO,  LENGTH,  NTOT,  NMODE,  NMODEB,  IPRINT,  KR.  DT,  RDT , 
WORK,  RWORK) 

Functional  Description 

Subroutine  MATMUL  computes  the  z-  and  y-forces  on  all  slender  body 
elements  and  for  all  modes  using  the  output  matrix  of  subroutine  BFSMAT, 
the  solution  matrix  obtained  by  subroutine  SOLVIT  and  the  C  aA,  C  aA 

pz  py 

columns  obtained  by  subroutine  MUZYC.  These  arrays  are  input  from 
utility  (scratch)  tapes.  As  the  matrix  multiplication  is  performed, 
the  resulting  body  forces  are  saved  in  column  order  on  two  utility 
tapes:  z-forces  on  tape  NEWBFZ,  y-forces  on  tape  NEWBFY. 


Input  Output  Variables 

MNEMONIC  DESCRIPTION 

NW  Tape  number  containing  the  solution  matrix 

(output  of  SOLVIT) 


NPSTAP 

NBFM 

NEWBFZ 

NEWBFY 

NWORK 

NTSBE 

NSBETO 


Tape  number  containing  the  CpzAA  and  Cp^AA  arrays 
Tape  number  containing  the  output  matrix  of 
subroutine  BFSMAT 

Tape  number  containing  the  output  z-forces 
Tape  number  containing  the  output  y-forces 
Dimension  of  the  work  array  WORK 
Length  of  the  CpzAA,  Cp^AA  arrays 
Total  number  of  slender  body  elements 
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2.16 


SUBROUTINE  HATS  (XA,  YA,  XO,  YQ,  LMBDA ,  X14C,  X34C,  Y,  H14C,  H34C, 

HP,  I HAT,  IFPRLL,  II 12,  NTH,  SCOEF,  THB) 

Functional  Description 

This  subroutine  computes  elements  of  the  h ^ /4Q- »  ^3/40”  anc*  ma^r“'x 
columns  for  lifting  surface  boxes  in  one  section  of  a  superpanel.  The 
matrix  SCOEF  calculated  in  subroutine  SPLN  can  be  used  to  determine  a 
function  ft  at  any  value  of  . 


|f(Ti)}  =  [sPLINEj^j  [Cn| 

=  jsPLINE^^j  [sCOEFj  I  PHI} 
»  [H]  {PHI} 


This  subroutine  outputs  [H]  at  the  1/4-  and  3/4-chord  points  of  each 
box  of  a  superpanel.  Also  output  is  the  derivative  of  [H]  at  the 

(nr 


3/4-chord  points  of  the  boxes.  Notice  that 


from 


SPLINE. 


(n) 


SPLINE.' 

J 


is  different 


The  former  is  for  the  nodal  points  i  =1,2- 


and  the  latter  is  for  the  interpolated  points  j  =  1,  2- 


Input  Output  Variables 


MNEMONIC 
XA,  YA 
XO,  YO 
LMBUA 
X14C 

X34C 


DESCRIPTION 

Inboard  edge  coordinates  of  elastic  beam 
Reference  coordinates  of  current  section 
A  =  sweep  angle  of  elastic  beam 
X-coordinate  array  of  the  1/4-chord  locations 
of  all  lifting  surface  boxes 
X-coordinate  array  of  the  3/4-chord  locations 
of  all  lifting  surface  boxes 
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MNEMONIC 


DESCRIPTION 


Y 

H14C 
H34C 
HP 
I  HAT 

IFPRLL 

NPAIR 

1112 

NTH 

SCOEF 

THB 


"Spanwise"  coordinate  array  of  all  lifting  surface 
boxes. 


Array  containing  a  column  of  the 


hl/4C 

h3/4C 

h' 


matrix 


Flag  to  select  h-,  a-»  or  0-calculations; 

IHAT  =  1  for  h,  2  for  a>  3  for  9 

Flag  to  indicate  sections  of  parallel  cutting  lines 

when  IFPRLL  =  1;  IFPRLL  =  0  otherwise 

The  number  of  II,  12  pairs  in  the  current  section 

of  a  superpanel . 

Array  of  the  'first-,  last-box  number  pairs' 

NTH  =  NODES+2  =  the  number  of  spline  coefficients 
for  the  current  superpanel 

Array  of  spline  coefficients  for  the  current  super¬ 
panels 

Array  containing  the  0B,  values  for  the  current 

J 

superpanel 


Calling  Subroutine  NEWH 
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mnemonic 

DESCRIPTION 

LENGTH 

Column  dimension  of  the  output  matrix  of 

subroutine  BFSMAT 

NTOT 

Number  of  unknowns 

NMODE 

Number  of  modes 

NMODEB 

Number  of  modes  on  bodies  only 

IPRINT 

Print  flag  for  body  element  forces 

IPRINT  =  0  for  no  print 

* 

KR 

Reduced  frequency 

DT 

Temporary  complex  work  array 

RDT 

Temporary  real  array  equivalent  to 

DT 

WORK 

Complex  work  array 

RWORK 

Temporary  real  array  equivalent  to 

WORK 

Calling  Subroutine  SDLM 
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.18  SUBROUTINE  MUZYC  (H,  DHDX,  K,  NTZY ,  NF BODY ,  NLBODY ,  NSBE,  KR, 


CBAR ,  AO,  AOP ,  XIS1,  XIS2,  AR,  NTPH ,  NTPDti,  UZY,  CPZY, 


UZYA ,  CPZYA ,  NR,  FMACH,  YB,  ZB,  CR,  NG 


Functional  Description 

Subroutine  MUZYC  calculates  the  axial  doublet  strengths  and  loading 
for  slender  bodies  both  in  symmetric  and  antisymmetric  conditions. 


Input  Output  Variables 


MNEMONIC 

H 

DHDX 

K 

NFBODY 

NLBODY 

NSBE 

KR 

CBAR 

AO 

AOP 

XIS1 

XIS2 

AR 

NTPH 

NTPDH 

UZY 


DESCRIPTION 

Array  containing  a  column  of  the  h-matrix 
Array  containing  a  column  of  the  dh/dx-matrix 
Flag  for  types  of  mode;  K  =  2  for  z-direction 
K  =  3  for  y-direction 
Number  of  first  body 
Number  of  last  body 

Array  of  the  number  of  slender  body  elements 
per  body 

Reduced  frequency 
Reference  chord 

Array  of  the  slender  body  element  half  widths  (radii) 
x-derivative  of  the  AO  array 
Leading  j 

/  edge  location  of  the  slender  body  elements 
Trailing  ) 

Array  of  the  cross-sectional  aspect  ratios  of  bodies 
Tape  number  (scratch  unit)  containing  the  h  matrices 
Tape  number  (scratch  unit)  containing  the 
dh/dx  matrices 

A  column  of  the  uz  or  Uy  array 
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MNEMONIC 


DESCRIPTION 

A  column  of  the  CpZAA  or  Cp^AA  array 
A  column  of  the  Uz  or  array  for  antisymmetric 


CPZYA 


A  column  of  the  CpzAA  or  CpyAA  array  for  anti¬ 
symmetric  gust 
Flag  for  gust  calculations 

Array  of  the  direction  cosines  for  gust  calculations 
Number  of  gusts 


FMACH,  YB,  ZB 
H,  DHDX 


See  Section  VI-3 

Arrays  of  the  h  and  dh/dx  matrices 


Calling  Subroutine  SB 
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2.19  SUBROUTINE  NEWH  (NTP5.  NTP6,  NTAPE,  MASTH4,  MASTH.  MASTDH,  NTOT, 

NVSM.  NASYM.  NROW,  MODES,  NBARAY,  NP,  NCARAY ,  YS,  ZS, 

Y,  X14C,  X34C .  H34C,  HP  NODE,  XI,  ETA,  ELXI,  ELYI , 

ELZI,  HCOL ,  PHI ,  THB,  IFP,  IPAIR,  1112,  SWQRIC,  WORK, 

NiiOETO,  NDOFTO,  PH  IT,  MDOF,  NODES) 

FUNCTIONAL  DESCRIPTION 

Subroutine  NEWH  is  the  working-main  of  the  h-,  a-,  e-calculation 
subroutines  for  lifting  surcace  panels  only.  It  reads  card  input  for 
all  superpanels  and  calls  subroutine  SPLN  for  the  calculation  of  the 
spline  coefficients;  it  also  calls  subroutine  HATS  which  calculates 
the  slopes  and  deflections  at  all  aerodynamic  elements  due  to  a  unit 
deflection  of  each  degree  of  freedom  at  each  node,  taken  one  at  a  time. 
That  is,  HATS  produces  deflection  and  slope  influence  coefficient 
matrices  [H],  [H^  [dH/dx].  These  influence  coefficient  matrices 
are  used  along  with  the  modal  matrix  [PHI],  in  NEWH,  to  produce  the 
slopes  and  deflections,  [h^Chy^J,  [dh/dx]  due  to  the  modes 

[h]  =  [H]  [PHI] 

^h1/4c^  =  £H1/4c^  [pHI^ 

[dh/dx]  =  [dH/dx]  [PHI] 

Also  in  NEWH  the  [PHI]  matrix  is  extracted  from  the  general  modal 
matrix  [PHIT],  The  [PHIT]  matrix  contains  up  to  eight  degrees  of 
freedom  per  lifting  surface  superpanel  node,  i.e.,  h,  »,  e,  f,  t,  4s 
6,  c  where  as  [PHI]  contains  only  the  up  to  three  values,  i.e.  h  a  e. 
Also  special  procedures  place  the  control  surface  or  tab  degrees  of 
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I 


freedom,  6  or  6,  into  the  a  degrees  of  freedom  for  a  control  surface 
superpanel . 

The  slopes  used  above  are  taken  at  the  3/4-chord  point  of  each  box. 

The  deflections  are  taken  at  the  same  point,  [h],  and  also  at  the  1/4- 
chord  point  [h1/4(J.  These  matrices  are  output  on  three  tapes;  [h] 
is  on  MASTH,  [dh/dxj  is  on  MASTDH  and  [h-)/4c]is  on  MASTH4  . 


Input  Output  Variables 


MNEMONIC 
NTP5 ,  NTP6 

NT  APE 

MASTH4 

MASTH 

MASTDH 

NTOT 

NSYM 

NASYM 

MOOES 

NBARAY 


DESCRIPTION 

Reference  number  of  the  system  input/output 
data  sets  (NTP5  =  5,  NTP6  =  6  on  IBM  370) 
Direct  access  unit  on  which  the  h,  »,  e 
columns  are  stored  (temporarily) 

Tape  number  containing  the  h^4c-columns 
Tape  number  containing  the  hg/^-columns 
Tape  number  containing  the  dh/dx-columns 
Total  number  of  lifting  surface  boxes 

Total  number  of  | SJESj^tri c ( 
modes 
MOOES 


nbi  = 


NSYM  +  NASYM 


|  _  1  nck  nsk,  i  =  1,  NP,  where 


NP  =  number  of  panels 


nc  =  number  of  chordwise  boxes  | 

*  In  panel 

nsk  =  number  of  spanwise  strips) 


NCARAY 


Array  of  the  number  of  chordwise  boxes  per  panel 
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MNEMONIC 
YS,  ZS 
Y 

X14C 

X34C 

H14C  f 

H34C  ■ 

HP 

NODE 

XI,  ETA 

ELXI  \ 

ELYI  ; 

ELZI  ) 

HCOL 

PHI 

THB 

IFP 


DESCRIPTION 

Arrays  of  the  y-  and  z-coordinates  per  strip 
Array  of  the  y-coordinates  of  all  boxes 
Array  of  the  1/4-chord  x-coordinates  of  all 
boxes 


Array  of  the  3/4-chord  x-coordi nates  of  all 


boxes 


Array  containing  a  column  of  the 


(hl/4c| 

lh3/4cf 


matri x 


Array  containing  a  column  of  the  dh/dx-(i .e.h ' ) 
matrix 

Array  containing  the  input  node-numbers  for 
any  one  superpanel 

S-,  n-arrays  (x-and  "spanwise"  coordinates)  of 
all  nodes  in  any  one  superpanel 


Arrays  •'n  /  of  all  nodes  for  the  case 

U ) 

Array  for  the  temporary  storage  of  any  one 
of  the  h  (or  h')  arrays 
Modal  matrix 

Array  containing  the  eB-column  (see  subroutine 
SPLN) 

Array  containing  the  "parallel-flags"  of  all 
sections  of  a  superpanel 
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MNEMONIC 


DESCRIPTION 


IPAIR 


1112 

SWORK 

WORK 

NODE TO 

N DO FTC 

PHIT 

KDOF 

MODES 


Array  containing  the  number  of  "first,  last" 
box-pairs  in  a  section  for  all  sections  of  a 
superpanel 

Array  containing  the  first,  last  box  numbers  for 
all  sections  and  all  superpanels 
Work  array 

An  NROW  by  NROW  work  array 
Total  number  of  modes 
Total  number  of  degrees  of  freedom 
The  general  modal  matrix 

Maximum  number  of  degrees  of  freedom  per  node 
Number  of  nodes  on  all  superpanels 


Calling  Subroutine  SDIM 

Called  Subroutines  ATAN3,  CTLS,  HATS,  SPLN,  ZEROUT 
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2.20  SUBROUTINE  0R6N(NTAPE,  MASTAP,  NSYM,  NASYM,  NODES.  MODES,  I ROW, 
COL,  PHI,  WORK) 


Functional  Description 

Subroutine  SPLINE  has  calculated  an  interpolation  matrix  for  unit 
deflections  at  nodal  points  taken  one  at  a  time,  i.e.,  an  interpolation 
influence  coefficient  matrix.  To  obtain  the  deflection  for  a  set  of 
modes,  [$],  a  matrix  multiplication  must  be  done. 

[H1/4]  »  [h1/4]  [*],  [fi3/4]  =  [h3/4]  [*],  [dH/dx]  =  [dh/dx]  M 


These  multiplications  must  be  done  for  both  symmetric  and  antisymmetric 
matrices  and  for  lifting  surfaces  and  bodies.  Subroutine  ORGN  performs 


these  matrix  multiplications. 
Input  Output  Variables 


MNEMONIC 

NTAPE 

MASTAP 

NSYM 

NASYM 

N0DES 

M0DES 

I  ROW 

COL 

PHI 

WORK 


DESCRIPTION 

Input  tape  containing  one  of  the  matrices  [H] 

Output  tape  containing  one  of  the  matrices  [H] 

Number  of  symmetric  modes 

Number  of  antisymmetric  modes 

Number  of  nodes 

Total  number  of  modes 

Maximum  number  of  row  elements  in  the  [H]  matrices 
Array  for  the  temporary  storage  of  one  column 
of  the  [H]  matrices 

One  column  of  the  modal  deflection  matrices  [PHI] 
Work  array  containing  one  of  the  [H]  matrices 


Call ing  Subroutine  DOTP 
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2.21  SUBROUTINE  PISTON  (NB.  NBOX.  NSBETO,  KR ,  CBAR,  FMACH ,  NTPH,  NTPDU, 
NTPH4,  MASTCP  MASTFZ,  MASTFY ,  NEWTPH,  NBARAY,  NCARAY, 

NSBEA,  YB,  ZB,  AO,  AOP,  XIS1,  XIS2,  CG,  SG,  EE,  YS,  ZS, 

H,  DHDX,  HP.,  XC,  JRUN,  IPRINT ,  DCP,  DCPS,  DCPA,  DFQ, 

DFQS,  DFQA,  I  POINT,  NR,  CR) 

Functional  Description 

This  subroutine  is  called  from  subroutine  SDLM  only  if  the  input  item 
NKP  is  different  from  zero.  In  this  case  subroutine  PISTON  is  called 
once  for  each  reduced  frequency  in  the  array  for  which  piston  theory 
calculations  are  specified.  It  computes  pressures  on  lifting  surface 
boxes,  and  body  z-,  y- forces,  for  all  modes  and  all  gusts.  The 
pressure  columns  are  saved  on  tape  MASTCP,  while  the  body  z-forces 
(y- forces)  are  saved  on  tape  MASTFZ  (MASTFY).  An  input  IPOINT,  tells 
whether  "average"  or  "point"  pressures  and  forces  are  calculated  for 
the  blast  onset  flows.  (See  Section  II. 4  for  details). 

The  equations  implemented  in  this  subroutine  are  as  follows.  Equation 
(61)  is  used  for  lifting  surface  elements  and  equation  (80)  for  body 
elements,  for  the  case  of  motion  dependent  forces.  Equation  (64)  and 
(70)  are  used  for  lifting  surface  elements  and  equation  (84)  is  used 
for  body  elements,  for  the  case  of  blast  gust  dependent  forces. 

Input  Output  Variables 

MNEMONIC  DESCRIPTION 

NB 

NBOX 


NSBETO 


Number  of  bodies 

Total  number  of  lifting  surface  boxes 
Total  number  of  slender  body  elements 


MNEMONIC 

KR 

CBAR 

FMACH 

NTPH 

NTPDH 

NTPH4 

MASTCP 

MASTFZ 

MASTFY 

NEWTPH 

H 

DHDX 

HP 

K 


JRUN 

IPRINT 

DCP,  DCPS 


DESCRIPTION 

Reduced  frequency 
Reference  chord 
Mach  number 

!h  arrays  \ 
dh/dx  arrays) 

hI/4c  arrays) 

Tape  number  on  which  the  pressures  on  lifting 
surface  boxes  are  saved 

Tape  number  on  which  the  body ^z- forces  are  sav?d 

ly-forces  are  saved 

Tape  number  on  which  the  h-j^  -arrays 
Array  containing  a  column  of  the  [h]  matrix 

Array  containing  a  column  of  the  [dh/dx]  matrix 

Array  containing  a  column  of  the  [h*]  matrix 

A  four- dimensional  array  containing  the  x- 

coordinates  of  the  four  corners  of  all  lifting 
surface  boxes. 

Sequence  number  of  the  current  frequency  for 
which  piston  theory  is  used  (1  <_  JRUN  £  NKP) 
Detail  print  flag;  IPRINT  >  1  prints  all 
pressures  and  body  element  forces 
Complex  array  containing  a  column  of  the 
pressures  (ACp)  on  lifting  surface  boxerf  due 
to  either  a  mode  or  a  symmetric  gust  column. 
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MNEMONIC 

DCPA 


DFQ,  DFQS, 


DFQA 


IPOINT 
NR,  NG 
CR 


DESCRIPTION 

Complex  array  containing  a  column  of  the 
pressures  (ACp)  on  lifting  surface  boxes  due 
to  an  antisymmetric  gust  column 
Complex  array  containing  a  column  of  either  the 
z-  or  the  y-body  forces  due  to  a  mode,  or  the 
z-body  forces  due  to  a  symmetric  gust  column. 
Complex  array  containing  a  column  of  the  z-body 
forces  due  to  an  antisymnetric  gust  column. 1 

Note  that  y-bosy  forces  due  to  gust  columns  are 
stored  (temporarily)  in  arrays  DCPS  and  DCPA 
Flag  to  select  point  pressures  when  IPOINT  f  Q 
Number  of  gust  columns 

The  3  X  20  array  containing  the  gust  direction 
cosines 


The  following  variables  are  described  under  Basic  Variables  in  Section 
VI 1-3 : 

NBARAY,  NCARAY,  NSBEA,  YB,  ZB,  AO,  AOP,  XIS1 ,  XIS2,  CG,  SG,  EE,  YS, 

zs. 


Calling  Subroutine  SDLM 
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2.22 


SUBROUTINE  ROWDYZ(NFB,  NLB,  ROW.  NTZYS,  D,  DA,  DX,  DY,  DZ,  IRB, 

BETA,  IDZDY ,  NTAPE,  NTAP2 ,  SGR,  CGR,  IPRNT ,  NT12, 

NSBE,  AO,  YB,  ZB,  AR,  XIS1,  XIS2) 

Functional  Description 

This  subroutine  performs  the  logic  required  to  set  up  the  argument  list 
for  subroutine  DZY  for  the  purpose  of  calculating  a  row  of  the  D?  or  D^ 
matrix  both  for  symmetric  and  antisymmetric  conditions  . 

Input  Output  Variables 

MNEMONIC  DESCRIPTION 


NFB 

First  body  number 

NLB 

Last  body  number 

ROW 

Row  number  of  the  Dz  or  D^  matrix  being  calculated 

NTZYS 

Number  of  columns  to  be  calculated 

D 

A  row  of  the  Dz 

,  or  Dy  matrix  for  symmetry 

DA 

A  row  of  the  Dz 

,  or  Dy  matrix  for  antisymmetry 

DX 

x -coordinate 

DY 

y -coordinate 

of  the  receiving  point 

DZ 

z -coordinate  . 

BETA 

6  =  /  1  -  M2  , 

M  =  Mach  Number 

IDZDY 

Flag  required  for  subroutine  FLLD 

NTAPE 

Tape  number  (scratch  unit)  on  which  rows  of  the 

symmetric  matrix  D  are  saved 

NTAP2 

Tape  number  (scratch  unit)  on  which  rows  of  the 

antisymmetric  matrix  DA  are  saved 

SGR 

Sine  ) 

;  of  the  receiving  point  dihedral  angle 

CGR 

Cosine  J 

IPRNT 


Print  flag 
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The  following  variables  are  described  under  Basic  Variables  in  Section 
VI 1-3 : 

NT1 2 ,  NSBE,  AO,  YB,  ZB,  AR,  XIS1 ,  XIS2. 

Calling  Subroutine  DZYMAT 

Called  Subroutine  DZY 


i 
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2.23  SUBROUTINE  RWREC(IFLAG,  NTAPE,  A,  NCWORD,  NUMBR,  NEWTAP) 

Functional  Description 

This  subroutine  reads  and/or  writes  unformatted  complex  records  A  of 
dimension  NCWORD  on  tape  according  to  three  options  defined  by  the 
flag  I FLAG. 


Input  Output  Variables 


MNEMONIC 
I  FLAG 


NTAPE 

A 

NCWORD 

NUMBR 


NEWTAP 


DESCRIPTION 

Option  flag; 

IFLAG  =  0,  write  A  on  NTAPE, 

IFLAG  =  1,  read  A  from  NTAPE, 

IFLAG  =  2,  copy  A  from  NTAPE  onto 
NEWTAP,  'NUMBR'  of  times 
Input/output  unit  number  to  use 
Complex  array  to  be  read  and/or  written 
Length  of  array  A 

Number  of  arrays  to  copy  from  NTAPE  onto 
NEWTAP  when  IFLAG  =  2; 

NUMBR  =  0  otherwise 
Output  unit,  when  IFLAG  =  2; 

NEWTAP  =  0  otherwise 


Calling  Subroutines  BFSMAT,  GENF,  MATMUL,  SB,  SDLM 
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24  SUBROUTINE  SB  (NM,  NSYM,  NASVM,  NTP6,  NTPH,  NTPDH,  IPRINT,  WORK, 

OW,  CPZYA,  UZYA,  NBEA,  NSBE,  NC,  NS,  AO,  YB,  ZB, 

AR,  AOP ,  XIS1,  XIS2 ,  X,  CG,  SG,  VP.  ?P.  CR.  NGl 

Functional  Description 

Subroutine  SB  calculates  the  normalwash  at  lifting  surface  boxes  and 
interference  body  elements  due  to  slender  body  elements  both  for 
symmetric  and  antisymmetric  conditions. 

Input  Output  Variables 

MNEMONIC  DESCRIPTION 


NM 

Total  number  of  modes  for  bodies 

NSYM 

Number  of  symmetric  modes 

NASYM 

Number  of  antisymmetric  modes 

NTP6 

Data  set  number  of  the  system  output  data  set 

NTPH 

Tape  number  (scratch  unit)  containing  the  h 

matri ces 

NTPDH 

Tape  number  (scratch  unit)  containing  the  dh/dx 

matri ces 

IPRINT 

Print  flag 

DU 

Array  containing  a  column  of  the  incremental 

normalwash  matrix,  [AW] 

CPZYA 

Array  containing  a  column  of  the[Cpz  A  A],  or 
[Cpy  A  A]  matrices 

UZYA 

Array  containing  a  column  of  the  [M2]  or  [M  ] 

matri ces 

WORK 

Complex  work  array 
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MNEMONIC 


DESCRIPTION 


CR  Array  of  the  direction  cosines  for  gust 

calculations 

NG  Number  of  gusts 

The  following  variables  are  described  under  Basic  Variables  in 
Section  VII-3: 

NBEA,  NSBE,  NC,  NS,  AO,  YB,  ZB,  AR,  AOP,  XIS1,  XIS2,  X,  CG,  SG, 
YP,  ZP. 

Calling  Subroutine:  SDLM 

Called  Subroutines:  DUMULT,  DZYMAT,  MUZYC,  READD,  RWREC 


2.25  SUBROUTINE  SDLM  (ELXI,  ELY  I ,  ELXI,  PHIN,  PHIZ,  PHIY,  NODES,  NSYM, 


NASYM,  MFIX1,  MFIX2,  IFNEWH,  MDOF,  PHI,  NAS,  NASB,  NBARAY, 
NCARAY ,  NSARAY ,  ISSTR,  NSSTR,  I  FLA,  NBEA,  NT12,  NSBEA,  YB, 
ZB,  ARB,  AVR,  XLE,  XTE,  RIA,  TH1A,  AO,  AOP,  XIS1,  XIS2, 

CG,  CS,  EE,  SG,  YS,  ZS,  XIJ,  YIN,  ZIN,  COORD,  X,  XIC, 
DELX,  XLAM,  H,  DHDX,  DELA,  DT,  DTA,  RHS ,  XC) 

Functional  Description 

Subroutine  SDLM  is  essentially  the  main  program  of  the  entire  aero 
program  segment  of  VIBRA-6.  It  first  calls  the  subroutines  that 
generate  the  basic  geometry  arrays  and  deflections,  then  the  rest  of 
the  subroutines  are  called  in  an  overall  DO  loop  on  all  reduced 
frequencies. 

Input  Output  Variables 
MNEMONIC 
ELXI 
ELYI 
ELZI 
PHIN 
PHIZ 

PHIY 

NODES 
NSYM 
NASYM 
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DESCRIPTION 
x-coordinate  j 

y-coordinate  ,  array  of  the  modes 
z-coordinate  / 

Modal  deflection  matrix  for  all  nodes  on  panels 
Modal  deflection  matrix  for  all  nodes  on  bodies 
in  z-di recti  on 

Modal  deflection  matrix  for  all  nodes  on  bodies 

in  y-direction 

Maximum  number  of  nodes 

Number  of  symmetric  modes 

Number  of  antisymmetric  modes 


l 


MNEMONIC 


DESCRIPTION 


MFIX1 

MFIX2 

IFNEWH 


MDOF 

PHI 

H 

DHDX 

DT 

DTA 

RHS 

xc 


Two  mode  numbers  for  which  aerodynamic 
parameters  are  calculated 
Flag  to  select  degrees  of  freedom;  IFNEWH  =  0 
for  h-  alone,  IFNEWH  f  0  for  use  of  available 
sectional  data 

The  number  of  degrees  of  freedom  In  the 
sectional  data  U  -  matrix) 

The  total  modal  deflection  matrix,  [$>] 

Real  array  for  a  column  of  the  h  matrix 
Real  array  for  a  column  of  the  dh/dx  matrix 
Complex  array  for  a  row  of  the  DT  matrix 
Complex  array  for  a  row  of  the  DTA  matrix 
Complex  array  for  temporary  use 
X-coordi nates  of  lifting  surface  box  corners 


The  following  variables  are  described  under  Basic  Variables  in 
Section  VII-3: 

NAS,  NASB,  NBARAY,  NCARAY ,  NSARAY ,  ISSTR,  NSSTR,  I  FLA,  NT12,  NSBEA, 
YB,  ZB,  ARB,  AVR,  XLE,  XTE,  RIA,  TH1A,  AO,  AOP,  XIS1 ,  XIS2,  CG,  CS, 
EE,  SG,  YS,  ZS,  XIJ,  YIN,  ZIN,  COORD,  X,  XIC,  DELX,  XLAM,  DELA. 


Calling  Subroutine  CSDLM 

Called  Subroutines  GEOM,  SPLINE,  DOTP,  GEND,  BFSMAT ,  SB,  WANDWT, 

GUST,  RHSIDE ,  SOLVIT,  MATMUL,  GENF,  AERO,  RWREC , 
HEADNG,  NEWH ,  PISTON,  ZEROUT 
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MNEMONIC 


DESCRIPTION 


11 

12 
NXQ 
LMAX 
X1L 
X2L 

X3L 

YL 

XP 

YP 


Index  of  first  element  in  the  coordinate  arrays 
Index  of  last  element  in  the  coordinate  arrays 
Number  of  nodes 

Maximum  dimension  of  the  coordinate  arrays 
3/4-cliord  -j 

V  x-coordinates  of  panel  boxes 
1/4-chord  J 

x-coordinates  of  body  element  endpoints 

y-coordinates  of  panel  boxes 

The  extracted  x-coordinate  array  for  the  current 

superbody/superpanel 

The  extracted  y-coordinate  array  for  the  current 
superpanel;  YP  =  0  for  superbodies 


H  One  column  of  the  [h]  matrix 

DHDX  One  column  of  the  [dh/dx]  matrix 

ISP  Not  used 

NSBP  Not  used 

i 

NBCUM  nbai  =  £  (NSBEb  +  1)  for  all  bodies,  i  =  1 ,  NB 

b=l 

i 

NPCUM  npai  =  [  (ncp)(nsp)  for  all  panels,  i  =  1 ,  NP 

p=l 

where  ncp  =  the  number  of  chordwise  boxes  in  panel  p 
nsp  =  the  number  of  spanwise  strips  in  panel  p 
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2.26  SUBROUTINE  SORTflFT,  IF2,  NSUP,  INSUP,  IT,  12,  NX 


LMAX,  X1L ,  X2L ,  X3L ,  YL,  XP,  YP,  H,  DHDX , 

ISP,  NSBP,  NBCUM,  NPCUM,  NTPH,  NTPDH, 

NTPSH,  KPRINT,  NTP6,  KCUM1  ,  KCUM2 ,  KCUM3 , 

NTOT ,  I_SORT ,  NOSE.  SlIPRH .  SUPRDR) 

Functional  Description 

Subroutine  SORT  has  a  dual  role:  it  extracts  the  coordinates  of  all 
points  in  the  current  superbody/superpanel  when  it  is  called  from 
subroutine  SPLINE  with  the  flag  ISORT  set  to  1,  while  on  the  second 
call,  when  IS0RT=2,  it  expands  the  hand  dh/dx  columns,  generated  by 
subroutine  SPL2  for  one  superbody/superpanel,  into  full  h  and  dh/dx 
columns  and  saves  these  on  the  output  tape  NTPSH. 


Input  Output  Variables 

MNEMONIC 

IF1 

IF2 


NSUP 

INSUP 


DESCRIPTION 

Flag  to  select  either  surface  spline 
(IF1  =  0)  or  linear  spline  (IF1  =  1) 

Flag  to  choose  1 /4-chord  x-coordinates  for 
superpanels,  when  IF2  =  1. 

IF2  =  0  means  that  3/4-chord  x-coordinates 
are  used  for  superpanels,  and  body  element 
endpoints  for  superbodies. 

Number  of  bodies/panels  in  superbody/superpanel 
Array  of  bodies/panels  in  superbody/superpanel. 
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t 


MNEMONIC 

NTPH 

NTPDH 

NTPSH 

KPRINT 

NTP6 

KCUM1 

KCUM2 

KCUM3 

NTOT 

I  SORT 


SUPRH 

SUPROH 

NODE 

Calling  Subroutine 


DESCRIPTION 

Tape  number  containing  the  h  columns  for  a 
superbody/superpanel 

Tape  number  containing  the  dh/dx  columns  for  a 
superbody/superpanel 

Tape  number  containing  the  expanded  (full)  h  and 
dh/dx  arrays 

KPRINT  =  1  prints  the  h  and  dh/dx  arrays; 

KPRINT  =  0  no  print 

The  data  set  number  of  the  system  output  data  set 
Column  number  of  the  h-arrays 
Column  number  of  the  dh/dx-arrays 

Column  number  of  the  h-arrays  for  1 /4-chord  x-coordinates 
Length  (dimension)  of  the  expanded  h  and  dh/dx  arrays 
Flag  that  selects  one  of  the  two  segments  in  subroutine 
SORT; 

I SORT  =  1  branches  to  segment  1 , 

the  coordinate-extractor, 

ISORT  =  2  branches  to  segment  2, 

the  h-  and  dh/dx-expander 
Array  containing  a  column  of  the  [h]  matrices 
Array  containing  a  column  of  the  [dh/dx]  matrices 
Array  of  node  numbers 

SPLINE 


I 
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2.27  SUBROUTINE  SPLINE  (MASTH,  MASTDH,  MASTH4.  NHODEB,  NMODE,  ELXI. 

ELYI ,  ELZI ,  WORK  NSBEA,  NBARAY,  NCARAY,  NSARAY ,  XIS1, 

XIS2,  X,  CG.  SG,  YS,  ZS,  XIC,  PHI,  PHIZ,  PHIY,  NDOF, 

N0DES,  NM0DES,  NDOFB  ,  IFNEWH,  XIJ,  NBCUM,  INSUP,  NODE, 

XIQ,  ETAQ,  CSG,  X3L,  YL,  XP  ,  YP ,  H,  DHDX) 

Functional  Description 

The  four  subroutines  SPLINE,  SORT,  SPLl,  SPL2,  constitute  a  surface 
spline  interpolation  procedure  based  on  the  method  of  Harder  and 
Desmarais  (Reference  9).  The  equations  used  in  this  procedure  are 
documented  in  Section  1 1. 3.  A  set  of  interpolation  functions  are  used 
whose  coefficients  are  found  in  SPLl.  The  slopes,  dh/dx,  and  deflections, 
h,  are  then  determined  at  specified  points  (XP,  YP)  in  SPL2  using  these 
functions  and  coefficients.  Subroutine  SPLINE  organizes  the  arguments 
for  the  three  called  subroutines  SORT,  SPLl  and  SPL2.  It  also  saves 
the  final  h  and  dh/dx  matrices  on  tapes  for  all  body  nodes  (if  any), 
then  for  all  panel  nodes. 

Input  Output  Variables 

MNEMONIC  DESCRIPTION 


MASTH 

Tape  number  containing  the  [h]  matrices  for 

bodies  and  panels 

MASTDH 

Tape  number  containing  the  [dh/dx]  matrices  for 

bodies  and  panels 

MASTH4 

Tape  number  containing  the  [hj^c]  matrix  for 

panels 

NMODEB 

Number  of  columns  in  the  [h]  matrix  for  bodies 

NMODE 

Number  of  columns  in  the  [h]  Matrix  for  panels 
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MNEMONIC 


DESCRIPTION 


ELXI 

x-coordinate  v 

ELYI 

y-coordinate  '  array  of  the  nodes 

ELZI 

z-coordinate  / 

WORK 

Work  array 

PHI 

The  total  modal  deflection  matrix,  [*] 

PHIZ 

Modal  deflection  matrix  for  all  nodes  on  bodies 

in  z-di recti  on 

PHIY 

Modal  deflection  matrix  for  all  nodes  on  bodies 

in  y-di recti on 

MDOF 

The  number  of  degrees  of  freedom  in  the  ^  matrix 

NODES 

The  total  number  of  nodes 

NMODES 

The  total  number  of  modes 

NDOFB 

The  number  of  degrees  of  freedom  for  body  bays 

IFNEWH 

Flag  to  select  degrees  of  freedom;  see 

subroutine  SDLM 

1 

nb,  =  l  (nsbe  >1 ) ,  i  =  1 ,  NB 

1  j  =  T  3 

where  nsbe.  is  the  number  of  slender  body 

J 

elements  in  body  j 

NBCUM 

INSUP 

Array  of  the  body-nuneers/panel -numbers  in  the 

current  superbody/superpanel 

NODE 

Array  of  node  numbers  in  the  current  superbody/ 

superpanel 

X I Q 

x-coordinate  array  of  the  nodal  points  in  the 

current  superbody/superpanel 

hiCMONIC 

ETAO 

CSG 

X3L 

YL 

TP,  YP 

H 

DHDX 


DESCRIPTION 

y-coordlnate  array  of  the  nodal  points  In  the 
current  superbody/superpanel 
Cosines  of  the  modified  dihedral  angles  of 
stri ps 

Array  of  the  slender  body  element  end  points 
In  the  current  superbody 
y-coordlnate  array  of  all  elements  in  the 
current  superbody 

x,  y-coordlnates  of  all  lifting  surface  boxes 
In  the  current  superpanel 
A  column  of  the  h  matrix 
A  column  of  the  dh/dx  matrix 


The  following  variables  are  described  under  Basic  Variables  in 
Section  VII-3: 

NSBEA,  NBARAY ,  NCARAY,  NSARAY,  XI31,  XIS2,  X,  CG,  SG,  YS,  ZS,  XIC,  XIJ. 


Calling  Subroutine  SDLM 

Called  Subroutines  SORT,  SPL1  and  SPL2 
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2.28  SUBROUTINE  SPLN  (SPLINE,  SCOEF.  THB,  XI,  ETA,  XA,  YA,  XB.  YB. 
NODES,  NDIM,  IPRINT,  NTP6) 


Functional  Description 

This  subroutine  calculates  the  coefficients  to  be  used  in  HATS  for 
interpolation  of  h,  a  and  e  along  the  elastic  axis  (along  t).  Influence 
coefficient  type  of  modes  are  used.  That  is  coefficients  are  calculated 
for  unit  deflections  of  the  variable  to  be  interpolated,  one  node  at  a 
time.  The  unit  deflection  matrix  [T]  is  given  by  equation  (89)  and 
includes  a  linear  extrapolation  to  the  elastic  axis  end  point.  The  one 
dimensional  version  of  the  Harder  and  Desman's  spline  curve  used  is: 

N  ?  ? 

f  (t)  =  Cj  +  C2t  +  £  Cn  (r-Tn)  *n  (r-Tn) 

n=3 

The  coefficients  Cn  can  be  determined  if  ff-r)  is  Krown  at  N  ♦  2  points, 
i  .e. ,  »  1  *  1 »  2  N  +  2 


or,  |f(T.)}  =  [SPLINE  ^]  |CnJ 
But,  ]f(  =  [T]  |PHI| 

|CnJ=  [SPLINE]-1  [T]  jPHl( 

[SCOEF]  jPHl| 


The  subroutine  outputs  [SCOEF]  to  be  used  later  in  HATS. 


Input  Output  Variables 

MNEMONIC 

SPLINE 


SCjJEF 


XI,  ETA 


DESCRIPTION 

The  NDIM  by  NDIM  matrix  [SPLINE] 

The  NDIM  by  N0DES  matrix  of  spline  coefficients 
Array  containing  the  ®  values  *or  a  superpanel 
The  n-  arrays  for  a  superpanel 


MNEMONIC 

DESCRIPTION 

XA,  YA 

Inboard  i 

>  coordinates  of  elastic  beam 

XB,  YB 

Outboard  ) 

N(*DES 

Number  of  nodes  in  the  current  superpanel 

NDIM 

NDIM  =  N0DES  +  4 

IPRINT 

Print  flag;  IPRINT  =  2  prints  the  [SPLINE]  and 

[SC0EF]  matrices 

NTP6 

Reference  number  of  the  system  output  data  set 

Calling  Subroutine 

NEWH 

Called  Subroutine 

MIS1 
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SUBROUTINE  SPLHIF1,  NXQ,  LMAX,  NTP1  ,  NTP6,  XIQ,  ETAQ,  RHS,  XKD, 
IXCON ,  IYCON,  NMAX) 


Functional  Description 

Equations  36  and  39  are  solved  in  this  subroutine.  Consider  these 
equations  toqether: 


where 


[KD] 


This  set  of  equations  is  solved  using  unit  solutions: 


f  0 

!  A 


{A} 


f! 

I  F  I 


[CMAT]  {A} 


Placing  these  in  the  first  equation  gives: 


t: 


=  [KD]  [CMAT] 


This  subroutine  solves  for  CMAT.  Using  the  relation  between 
{a}  shows  that 


and 


[CMAT]  = 


See  equations  41b  and  42c. 


Input  Output  Variables 

IF1  Flag  to  select  either  surface  spline  (IF!  =  0) 

or  linear  spline  (IF1  =  1) 

NXQ  Number  of  nodes  in  the  current  superbody/superpanel 
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MNEMONIC 

DESCRIPTION 

LMAX 

Not  used 

NTP1 

Utility  (scratch)  tape  number 

NTP6 

The  data  set  number  of  the  system  output  data 

set 

XIQ 

x-coordinates  ) 

]  of  the  nodal  points 

ETAQ 

y-coordi nates  ) 

RHS 

Matrix  of  a  set  of  unit  deflections  for  which 

matrix  equation  is  solved  by  subroutine  MIS1 

the 

XKD 

Influence  coefficient  matrix 

IXCON 

I 

| 

|  Constant  flags  for  superbodies; 

|  IXCON  =  1,  x=constant;  IXC0N=0  otherwise; 

IYCON 

i 

|  IYCON  =  1,  y=constant;  IYC0N=0  otherwise 

NMAX 

The  row-  and  column-dimensions  of  the  square  matri 

RHS  and  XKD 

Calling  Subroutine  SPLINE 

Called  Subroutine  MIS1 
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2.30  SUBROUTINE  SPL2(IF1,  IF2,  NXQ,  LMAX,  NTP2,  NTP3,  NTP6,  KPRINT,  NTPH 


UTPUH ,  12.  XP.  YP,  XIQ,  ETAQ,  XKF,  DKFDX,  CMAT,  SUM,  NMAXj 


Functional  Description 

This  subroutine  finds  {h}  and  {dh/dx}  at  the  points  to  be  interpolated 
knowing  the  coefficients  [CMAT]  determined  in  subroutine  SPL1 .  This 
subroutine  calculates  [Ijk3  and  [3ljk/3Xj]  (see  Eqs.  (47)  and  (48), 
as  follows: 

[Ijk]  =  [KF]  [CMAT] 


where 


[DKFDX]  [CMAT] 


[KF]  »  [P.  Rjl] 

[OKFDX]  "  [SPj/SX  aRjt/sx] 


Input  Output  Variables 


MNEMONIC 

IF1 

IF2 


LMAX 

NTP2 

NTP3 


DESCRIPTION 

Flag  to  select  either  surface  spline  (IF1  =  0) 
or  linear  spline  (IF!  =  1 ) 

IF2  =  1  bypasses  the  dh/dx  calculations  when 
the  1 /4-chord  x-coordinates  are  used  for  panel 
boxes ; 

IF2  =  0  otherwise 

Maximum  dimension  of  the  aero  coordinate  arrays 

XP,  YP  for  the  current  superpanel 

Utility  (scratch)  tape  number 

Utility  (scratch)  tape  number 

The  data  set  number  of  the  system  output  data  set 


NTP6 


V 


MNEMONIC  DESCRIPTION 

KPRINT  Print  flag; 

KPRINT  =  1,  print  the  h  and  dh/dx  columns, 
KPRINT  =  0,  no  print 


NTPH 

Tape  number  containing  the  h  columns 

NTPDH 

Tape  number  containing  the  dh/dx  columns 

12 

Index  of  the  last  element  in  the  XP,  YP  arrays 

XP 

x- coordinate  1 

array  of  the  aero  points  in  the 

YP 

y-coordinate  ; 

current  superpanel /superbody 

X  IQ 

x- coordinate  ) 

array  of  the  nodal  points  in  the 

ETAQ 

y-coordinate  ) 

current  superpanel/superbody 

XKF 

The  KF-matrix 

DKFDX 

The  DKDF-matrix 

CMAT 

Matrix  of  the  set  of  coefficients  calculated  in 

subroutine  SPll 

SUM 

Result  of  the  matrix  multiplication. 

i .e. ,  one  of  the 

matrices  [h]  or  [dh/dx] 

NMAX 

Column  dimension 

of  the  matrix  SUM; 

Row  and  column  dimension  of  the  matrix  CMAT 

Calling  Subroutine  SPLINE 
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2-31  SUBROUTINE  SUBB(KB,  KS,  I,  J,  JZ,  J B,  LB,  L$,  NDY ,  NYFL,  FLND, 

FLNE,  PI,  EPS,  SGR,  CGR,  SGS,  CGS,  AR,  SL.  CL,  TL.  FL, 

BETA,  SUM,  YB,  ZB,  RIA ,  X,  YS,  ZS,  DELX) 

Functional  Description 

Subroutine  SUBB  computes  the  downwash  factor  matrix  elements  for  all 
receiving  points,  and  all  sending  points  on  interference  bodies,  one 
element  at  a  time,  both  for  symmetric  conditions  and  for  antisymmetric 
conditions.  The  results  of  the  calculations  are  returned  to  the  calling 
routine  in  the  array  SUM  via  the  argument  list. 


Input  Output  Variables 

MNEMONIC 

KB 

KS 

I 

J 

JZ 

JB 

LB 

LS 

NDY 

NYFL 

FIND 

FLNE 

PI 


DESCRIPTION 

Index  of  receiving  body.  0  when  receiving 
point  is  on  panel 

Strip  number  in  which  receiving  point  lies 
Receiving  point  index  (row  number  of  DT  matrix) 
Sending  point  index  (coiumn  number  of  DT  matrix) 
Not  used 

Sending  body  element  number 
Sending  body  number 

Index  of  sending  element  y-  and  z-coordinates 
Body  flag;  NDY  =  1  for  sending  points  in  y- 
oriented  bodies,  NDY  =  0  otherwise 
NYFL  =  0  for  DYP  and  DYZ  elements, 

NYFL  =  1  for  DYY  elements 
Not  used 

Ground  effect  flag 
"  *  3.14159 
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MNEMONIC 


DESCRIPTION 


EPS 

SGR 

CGR 

SGS 

CGS 

AR 

SL 

CL 

TL 

FL 

BETA 

SUM 


E  =  O.OOOOI 
Sine  j 

/  of  receiving  point  dihedral  angle 
Cosine  j 

Sine  | 

Cosine  j  of  sending  point  dihedral  angle 

Cross-sectional  aspect  ratio  of  receiving  body 

sinA  =  0 

cosa  =  1 

tanA  =  0 

Reference  chord 

B  =  \'1  -  M^ ,  where  M  =  Mach  Number 

One  element  of  the  DT  matrix,  and 
one  element  of  the  DTA  matrix,  for 
sending  points  on  an  interference  body 


The  following  variables  are  described  under  Basic  Variables  in  Section  VI 1-3 
VB,  ZB,  RIA,  X,  YS,  ZS,  DELX . 


Calling  Subroutines  DZPY,  DYPZ 

Called  Subroutine  DZY 


Subroutine  SUBP  computes  the  downwash  factor  elements,  DP^ ,  for  symmetric 
conditions,  and  DPA^j  for  antisymmetric  conditions,  for  all  receiving  points 
i  on  panels  and  interference  body  elements,  and  sending  points  j  on  panels, 
one  element  at  a  time.  The  results  of  the  calculations  are  returned  to 
the  calling  routine  in  the  array  SUM  via  the  argument  list. 


Input  Output  Variables 


MNEMONIC 

I 

L 

LS 

J 

10 

IR 

NBXS 


DESCRIPTION 


NCPNB 


'JiP 

,P 


in  which  the  sending  point  lies 


Receiving  point  index 
Panel  number 
Strip  number 
Column  number  of  DP  matrix 

Chordwise  box  number  that  corresponds  to  sending 
point  j 

Sending  point  index 

The  number  of  boxes  on  the  panel  in  which  the  sending 
point  lies  plus  the  total  number  of  boxes  of  the 
preceding  panels 

The  number  of  boxes  in  the  first  strip  of  the  panel 
in  which  the  sending  point  lies  plus  the  total 
number  of  boxes  of  the  preceding  panels 
sine  | 


of  the  dihedral  angle  of  the  receiving  strip 


cosine 


MNEMONIC 


DESCRIPTION 


YREC 

y-coordinate  j 

/  of  receiving  points 
z-coordinate  ’ 

ZREC 

SUM 

One  element  of  the  DP  matrix  and  one  element 

of  the  DPA  matrix 

SCALER 

Nondimensional  image  radius 

The  following  variables  are  described  under  Basic  Variables  in 
Section  VI I -3 : 

NCARAY ,  YB,  ZB,  ARB,  AVR,  XLE,  XTE,  X,  CG,  EE,  SG,  YS,  ZS,  XIC, 
DELX,  XLAM  : 

Calling  Subroutines  DPPS,  DPZY 

Called  Subroutines  SNPDF,  INCRO,  SUBI 
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2.33  SUBROUTINE  WANDWT (IPRINT .  LINES.  NOUT,  NTQT,  NBQX,  NB,  NDW, 

NTAPH,  NTPDH,  NWT,  NWTA,  REFC,  KR,  H,  DHDX,  W,  DW,  COL) 

Functional  Description 

This  subroutine  calculates  the  complex  downwash  boundary  conditions,  W, 
on  the  lifting  surfaces  due  to  lifting  surface  motion,  and  combines  these 
with  the  normalwash,  aW,  due  to  the  slender  bodies.  Both  symmetric  and  anti¬ 
symmetric  motions  are  considered.  The  normalwash  matrix  for  symmetry,  WT, 
is  saved  on  tape  NWT,  the  one  for  antisymmetry,  WTA,  is  saved  on  tape  NWTA. 


Input  Output  Variables 


MNEMONIC 

IPRINT 

LINES 

NOUT 

NTOT 

NBOX 

NB 

NDW 

NTAPH 

NTPDH 

NWT 

NWTA 

REFC 

KR 


DESCRIPTION 

Print  flag 

Maximum  number  of  lines  per  page 

Data  set  number  of  the  system  output  data  set 

Number  of  unknowns 

Number  of  lifting  surface  boxes 

Number  of  bodies 

Tape  number  containing  the  AW  matrix 

Tape  number  containing  the  h  matrix 

Tape  number  containing  the  dh/dx  matrix 

Tape  number  on  which  the  WT  matrix  (normalwash 

columns  for  symmetry)  is  written 

Tape  number  on  which  the  WTA  matrix  (normalwash 

columns  for  antisymmetry)  is  written 

Reference  chord 

Reduced  frequency 
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MNEMONIC 


DESCRIPTION 


H 

A  column  of  the  h  matrix 

DHDX 

A  column  of  the  dh/dx  matrix 

W 

A  column  of  the  W  matrix 

DW 

A  column  of  the  AW  matrix 

COL 

Complex  array  for  temporary  storage 

Calling  Subroutine 

SDLM 

2.34  SUBROUTINE  ZEROUT(WORK,  LENGTH,  LOOP,  ITAPE) 

Functional  Description 
This  subroutine  initializes  the  complex  array  WORK  of  dimension  LENGTH 
to  zeroes.  Also,  this  array  of  zeroes  in  WORK  may  be  written  on  tape 
ITAPE  'LOOP'  times  whenever  the  arguments  LOOP  and  ITAPE  have  non-zero 
values, 

Input  Output  Variables 

MNEMONIC  DESCRIPTION 


WORK 

Complex  work  array 

LENGTH 

Length  of  array  WORK 

LOOP 

The  number  of  zero  arrays  to  be  written 

on  the  output  unit  ITAPE 

ITAPE 

Output  unit  number 

Calling  Subroutine 

GENF,  SB,  SDLM 
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3. 


BASIC  VARIABLES 


The  following  variables  are  used  throughout  the  aerodynamics  module 
of  program  VIBRA-6. 


MNEMONIC 

NB 

ND 

NE 

NP 

NBY 

NBZ 

NTO 

NTP 

NTY 

NTZ 

NTYS 

NTZS 

MAXGR 

MAXSTR 

NSBETO 

NSTRIP 

KR 

XM 

RE  FA 
REFC 


DESCRIPTION 

Number  of  bodies 
ND  =  0 

e ,  ground  effect  flag 

Number  of  panels 

Number  of  y-oriented  bodies 

Number  of  z-oriented  bodies 

NTO  =  NTP  +  NTZ  +  NTY  (see  below) 

Total  number  of  lifting  surface  boxes 
Number  of  y-oriented  interference  body  elements 
Number  of  z-oriented  interference  body  elements 
Number  of  y-oriented  slender  body  elements 
Number  of  z-oriented  slender  body  elements 
Number  of  components  (groups  of  panels) 

Number  of  superstrips  on  all  panels 
Total  number  of  slender  body  elements 
Number  of  lifting  surface  strips 
Reduced  frequency 
Moment  axis 

Reference  area;  usually  area  of  both  wings 
Reference  chord 
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MNEMONIC 


DESCRIPTION 


REFS 

FMACH 

LINES 

NAS 

NASB 

NBARAY 

NCARAY ,  NC 

NSARAY ,  NS 

ISSTR 
NSSTR 
I  FLA 

NBEA 

NT1 2 

NSBEA,  NSBE 

YB 

ZB 

ARB 

AVR 

XLE 


Reference  semi  span 
Mach  Number 

Maximum  number  of  print  lines  per  page 

Number  of  associated  bodies  per  panel 

Associated  body  numbers  for  all  panels 
i 

nba,  =  l  nc.  ns,,  i  =  1 ,  NP  (see  below) 

1  j=l  3  3 

nci  =  number  of  chordwise  boxes  per  panel, 
i  =  1,  NP 

nsi  =  number  of  spanwise  strips  per  panel, 
i  =  1,  NP 

Superstrip  number  of  each  strip 
Number  of  strips  per  superstrip 
Sequence  numbers  of  the  first  and  last  body 
element  per  body 

Number  of  interference  body  elements  per  body 
and  their  z-y  orientation  flags 
Number  of  u-points  (on  interference  body  elements) 
per  body 

Number  of  slender  body  elements  per  body 
y-coordinate  array  of  body  centerlines 
z-coordinate  array  of  body  centerlines 
Cross-sectional  aspect  ratio  per  body 
Average  characteristic  half-width  per  body 
Leading  edge  x-coordinate  per  body 
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MNEMONIC 


DESCRIPTION 


XTE 

RIA 

TH1A 

TH2A 

AO 

AOP 

XISI 

XIS2 

CG 

CS 

EE 

SG 

YS 

ZS 

XIJ 

YIN 

ZIN 

COORD 

X 

XIC 


Trailing  edge  x-coordinate  per  body 

Body  interference  element  radii 

Array  of  the  THl-values  (angular  orientation) 

for  all  bodies 

Zero  array 

Slender  body  element  half-widths  (y-direction) ; 
radii  for  circular  cross  sections 
x-derivative  of  the  AO  values  for  all  slender 
body  elements 

Leading  edge  x-coordinates  of  slender  body  elements 
Trailing  edge  x-coordinates  of  slender  body  elements 
Cosine  of  the  dihedral  angles  of  all  strips 
Chord-length  of  strips 
Half-width  of  strips 

Sine  of  the  dihedral  angles  of  all  strips 
y-coordinate  ) 

/  of  strip  centerlines  and  body  elements 
z-coordinate  ) 

x-coordinate-  of  leading  edge  of  strip  centerlines 
y-coordinate  of  inboard  edge  of  panel 
z-coordinate  of  inboard  edge  of  panel 
Spanwise  coordinate  of  strips 
3/4-chord  x-coordinates  of  all  lifting  surface 
bodies  and  interference  body  element  midpoints 
1/4-chord  x-coordinates  of  all  lifting  surface  boxes 
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MNEMONIC 


DESCRIPTION 


v„ 


DELX 

XLAM 

DELA 


Average  chord-lengths  of  all  lifting  surface 
boxes  and  interference  body  element  lengths 
Tangent  of  the  sweep  angle  of  the  1/4-chord 
line  (bound  vortex)  of  all  lifting  surface  boxes 
Lifting  surface  box  areas 
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APPENDIX  A 


THE  BLAST  ONSET  PRESSURE  FIELD 

The  analyses  of  References  13  and  14  show  that  the  loads  induced  by 
the  diffraction  of  a  blast  wave  about  a  thin  airfoil  are  equivalent  to 
those  obtained  when  only  the  travelling  gust  part  of  the  problem  is  con¬ 
sidered.  Quoting  Reference  14,  "It  is  the  purpose  of  this  note  to  show 
that  when  a  blast  wave  of  any  magnitude  is  represented  by  a  travelling 
gust  and  the  impulse  is  computed  using  apparent-mass  concepts,  the  re¬ 
sult  is  equivalent  to  accounting  for  diffraction  of  the  blast  wave".  The 
article  then  shows  that  the  apparent-mass  concept  (acoustic  theory)  lies 
within  10%  of  the  correct  value  for  blast  waves  ranging  in  strength  from 
acoustic  (weak)  to  infinitely  strong. 

The  above  conclusion  can  be  understood  more  clearly  if  basic  acousitc 
theory  is  considered.  First,  it  should  be  stated  that  the  gust  or  velocity 
field  and  the  overpressure  field  are  just  different  aspects  of  the  blast 
wave  and  are  not  independent.  The  velocity  field  can  be  obtained  from  the 
pressure  field  and  vice  versa.  Acoustic  problems  are  solved,  in  general, 
as  follows.  The  incoming  wave  creates  an  onset  flow  at  the  aircraft  surface. 
This  onset  flow  is  not  really  fundamentally  different  from  any  other  onset 
flow  field,  e.g.,  a  flow  at  an  angle  of  attack.  A  disturbance  flow  must  be 
generated  at  the  surface  so  that  the  sum  of  the  two  (i.e.,  the  onset  and 
disturbance  flow)  produces  streamlines  that  pass  around  the  body  and  thus 


13.  Drischler,  J.A.,  and  Piederich,  F.W.,  "Lifit  and  Moment  ReAponAeA  to 
Pe.ne.tsia.tion  o  ft  Sha/ip  -  Edged  Ttiavefing  Gua  tA ,  uxitk  Application  to  Pene.tA.a- 
tion  ofj  Weak  BlaAt  WaveA",  NACA  TN  3956,  May  1957. 

14.  Goodman,  T.R.,  and  Sargent,T.P. ,  "Some  Evidence  JuAtifiying  the  Reps te- 
ientation  oft  BlaAtA  a*  Traveling  GuAtA” ,  J.Aero.  Sci.,  Vol .  26,  No.  9, 

p.  608,  1959. 
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satisfy  the  body  boundary  condition  of  no  flow  through  the  body  surface.  The 
final  pressure  field  is  then  obtained  by  adding  the  pressure  fields  from  the 
onset  and  disturbance  flows.  For  a  lifting  surface  the  disturbance  pressure 
takes  a  jump  across  the  surface  at  all  times,  whereas  the  onset  flow  does  not. 
The  difference  in  pressure  across  the  surface  of  the  onset  flow  depends  on 
the  thickness  of  the  surface  and  the  gradient  of  the  pressure  field. 

The  conclusion  reached  in  Reference  14  is  that  the  onset  pressure  field 
does  not  affect  the  surface  to  any  noticeable  degree  and  that  the  distur¬ 
bance  pressure  field  (called  added-mass  effect)  contains  all  of  the  lifting 
effect.  It  is  postulated  that  the  reason  for  this  conclusion  depends  on 
the  fact  that  a  "thin"  airfoil  was  the  surface  considered.  To  be  completely 
accurate  the  lift  induced  by  the  onset  pressure  field  should  be  added, 
however  small  it  may  be. 
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Consider  Figure  37.  The  jump  in  pressure  across  the  surface  is  given 
in  terms  of  the  pressure  gradient  normal  to  the  surface  3P/an  and  the  sur¬ 
face  thickness  e  as  follows: 


L  iL 

q  an 


In  the  frequency  domain  the  overpressure  is  given  by: 


P  =  P0  exp  (-i2kRi/<T) 


where 

R  =  M  /(I  +  M  cosa) 

co  ao 

i  -  X  COSa  +  y  cose  +  Z  COSy 


and 


n  _ 

32 


cosy; 


3i 

ay 


cose 


and  cosa,  cose  and  cosy  are  the  direction  cosines  of  the  blast  orienta 
tion.  The  final  result  for  ACp  is  then, 

ACn  =  '  E  (cosy  cosy  -  ccse  siny)  exp  (-i2kRt/<T) 

P  Pc 


The  surface  loading  due  to  the  onset  pressure  field  is  proportional 
to  c,  the  surface  thickness,  and  is  therefor  usually  small  for  conven¬ 
tional  aircraft  wings.  Even  though  it  is  small,  however,  it  should 
be  retained  in  the  analysis. 


The  loading  across  a  fuselage  may  now  be  considered.  Let  aF/ax 
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V 


be  the  force  (normalized  by  dynamic  pressure)  on  a  fuselage  cross  section 
in  the  direction  of  motion  of  the  blast  wave.  This  force  is  the  result 
of  integrating  the  onset  pressure  field  around  a  circle  and  is  given 
by: 


Figure  A-Z  Body  Geometrical  Variables 


where,  in  Figure  3S,  z  is  in  the  direction  of  the  blast  and  y  is  parallel 
to  the  blast  wave  front  centered  on  the  body  axis.  The  quantities  of 
Interest  in  the  Integral  are: 

rf  •  F  =  (j  cose  +  k  sine)  •  F  =  sine 

P  =  P0  exp  (-i2kftt/c) 

T  *  X  COSa  +  Yc  COSB  +  Zc  COSY 

t  *  Z  COSa  +  T 

where 

COSy  »  Vcos^e  +  COS^Y 
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N- 


■  raw  'ni-v 


and  2  ■  a  sine.  The  quantities  Yc,  Zc  are  the  coordinates  of  the  fuse- 

sF 

lage  center  In  aircraft  coordinates.  Then  f-  becomes 

oX 

|£  -  ,  ^  exp  (-12MK/E)  I  2"  s,„9  de 

0 

or,  upon  noting  that  s1n(e  +  *}  =  -sine,  one  obtains: 

it.  !?exp  H2kRt/c)  [’  21  sine  de 

®x  H  'q  C 

The  integral  can  be  expressed  as  a  Bessel  function  since 

d  Jn  (*)  l  /  , 

— g-jj —  m  ~  *  I  sin  (X  sine  -  ne)  sine  de 

0 

Thus,  if  n  =  0  and  X  *  2kRa  cosy/c,  then  becomes: 

=  a  ~  exp  (-12kRT/rr)  2*1  fjjgM 

and  after  differentiating  the  Bessel  function,  the  final  result  Is: 

|£  B  _  a  _9.  2*1  J,  (2kRa  cosy/c)  exp  (-12kRl/c) 

dX  Cj  I 

The  force  per  unit  length  In  the  y  direction  (parallel  to  the  blast 
front)  is  zero  since  the  pressure  does  not  vary  in  that  direction.  Thus, 
the  force  per  unit  length  In  any  direction  can  be  obtained.  Specifically, 
the  running  loads  In  the  z-  and  y-dlrectlons  are: 

ll-iL  *  if  cosy 

3X  3X  COSy 

III .  if.  cose 

3X  3X  COSy 

Currently,  the  aerodynamics  module  calculates  only  loadings  due  to 
the  "disturbance  flow".  This  module  can  be  modified  to  produce  "onset 
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flow"  loadings  as  well.  The  "onset  flow”  loadings  would  be  passed  to  the 
frequency  response  and  unit  gust  load  modules  and  unit  solutions  generated 
In  the  same  fashion  as  presently  done  for  the  "disturbance  flow".  These 
solutions  would  then  pass  to  the  blast  and  time  response  modules  for 
multiplication  with  the  Fourier  transforms  of  the  gust  velocity  and  onset 
pressure  time  histories.  The  Inverse  transform  summation  of  the  two 
resulting  frequency  response  solutions  would  then  yield  the  time  history 
of  the  loads  and  responses  to  the  total  loading. 
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SYMBOLS 


A 

a 

AR 

b 

c 

c 


3x  3x 
f 


Reference  area  (full  area) 

Body  radius 

Aspect  ratio  of  body  cross-section  (=height/width) 
Ming  span 

Local  chord  length 
Reference  chord  length 

Rolling  moment  coefficient  (Rolling  Moment)/qAb 
Pitching  moment  coefficient  (Pitching  Moment)/qAc 
Lateral  force  coefficient  (Lateral  Force)/qA 
Vertical  force  coefficient  (Vertical  Force)/qA 
Doublet-Lattice  Method 
Strip  semi  width 

Running  load  (normalized  by  q)  on  bodies  in  y  and  z 
directions 

Frequency  in  HZ  or  cps.  f  =  7^ 

£7T 


G(t)  Gust  particle  velocity  profile  as  a  function  of 

distance  of  penetration,  s 

Gq  Fourier  transform  of  G(c) 

G  Symmetry  factor  for  lifting  surfaces  that  lie  in  the 

symmetry  plane 

g  Symmetry  factor  for  bodies  that  lie  in  the  symmetry 

plane,  also  gust  frequency  response 

h ,a ,6  Modal  deflections  at  a  modal  point,  h;  deflection 

normal  to  lifting  surface.  a,0;  angular  deflections 
pitch  and  bending  in  elastic  axis  coordinates 
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<» 


U.t 


alj^/aX 

o„0-),v„(P) 

kr 


M  or  M 

ao 

M 

n 

P 

PT 

q 

q 

s 

t 

U  or  U 

GO 

g 

w 

x.y.z 

xa ,ya,za 
a.S.Y 


Unit  vectors  in  the  x,  y  and  z  directions 

Spline  interpolation  matrix  whose  elements  relate  the 
deflection  of  point  j  due  to  the  deflection  at  point  k 

Slope  of  spline  interpolation  matrix 

Bessel  functions 

Reduced  frequency,  u>c/2U 

Distance  from  blast  to  point  on  aircraft  measured  in 
a  line  parallel  to  the  blast  direction  of  travel 

Free  Stream  Mach  number 

Added  mass  of  body  cross  section 

Unit  vector  normal  to  lifting  surface 

2kr  M^a/c 

Piston  Theory 

Generalized  modal  coordinates 

2 

Dynamic  pressure,  1  pU 

Semi  span,  b/2 
Time 

Frees tream  velocity,  or  aircraft  velocity 
Velocity  vector  of  gust  wave  front 
Unit  vector  in  direction  of 
Normalwash  on  lifting  surface 
Right-handed  coordinates  fixed  In  aircraft 
Right-handed  coordinates  fixed  at  blast  origin 

Direction  cosines  of  the  blast  wave  velocity  when  it 
hits  the  aircraft  in  aircraft  coordinates;  y  Is  also 
dihedral  of  panel 
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Difference  of  upper  and  lower  lifting  surface 
pressures 

Coefficient  form  of  AP 
Fraction  of  local  chord 
Angular  position  on  body  surface 
Elastic  axis  sweep 
Wave  length  ■  Z^^ 

Modal  matrix 

Distance  along  elastic  axis, 

also  t  =  2kr  R/c 

ax  kr  R  coso/c 
Frequency,  rad/sec 
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